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AFIT/GAE/ENY/95D-10 


Abstract 

Mach  2.9  boundary  layer  flow  (R./m  w  1.75x10^)  under  the  influence  of  mild  pressure 
gradients  is  studied  numerically.  Baldwin-Lomax  and  k  —  u  turbulence  models  are  incor¬ 
porated  into  a  cell-centered  finite  volume  flow  solver  and  the  results  are  compared  with 
hot  wire  anemometry  and  Laser  Doppler  Velocimetry  (LDV)  measurements  obtained  for 
the  same  geometries  in  the  AFIT  Mach  2.9  wind  tunnel.  Agreement  between  the  present 
simulations  obtained  with  the  A:  -  w  turbulence  model  and  experimental  velocity  profiles 
IS  excellent  in  all  test  sections.  Nondimensional  turbulent  shear  stress  predictions  closely 
match  experimental  data  in  the  flat  plate  and  adverse  pressure  gradient  sections  while 
slightly  over-predicting  this  quantity  in  the  favorable  pressure  gradient  region.  Favorable 
pressure  gradients  are  found  to  stabilize  the  flowfield,  resulting  in  increased  boundary 
layer  thickness  and  reduced  turbulent  and  wall  shear  stress  distributions.  Additionally, 
the  presence  of  a  favorable  pressure  gradient  is  found  to  diminish  the  effects  of  variations 
in  upstream  boundary  condition  specification.  Adverse  presure  gradients  are  found  to  de¬ 
stabilize  the  flowfield,  resulting  in  increases  in  the  turbulent  shear  stress,  turbulent  kinetic 
energy,  and  wall  shear  stress.  Upstream  effects  are  found  to  play  a  major  role  in  adverse 
pressure  gradient  flowfield  development.  Flowfield  features  are  predicted  more  accurately 
with  the  k  —  u  model  than  with  the  Baldwin-Lomax  model. 
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NUMERICAL  SIMULATION  OF  SUPERSONIC 


TURBULENT  BOUNDARY  LAYER  FLOW 
UNDER  THE  INFLUENCE  OF  MILD 
PRESSURE  GRADIENTS 

I.  Introduction 

1.1  Problem  Motivation 

The  vast  majority  of  flowfields  of  current  engineering  interest  are  turbulent  in  na¬ 
ture.  High  speed  turbulent  flows  are  of  particular  interest  to  the  Air  Force  as  many  of  the 
weapon  systems  currently  under  development  operate  in  the  supersonic  and  hypersonic 
regimes.  As  the  atmospheric  conditions  in  which  these  aerospace  vehicles  operate  become 
more  and  more  hostile,  experimental  determination  of  the  character  of  these  flowfields  be¬ 
comes  prohibitively  difficult  and  expensive.  In  the  current  atmosphere  of  reduced  defense 
spending,  the  emphasis  must  be  placed  on  working  smarter,  not  harder.  Numerical  simu¬ 
lation  of  these  flowfields  must  be  an  option  by  which  design  issues  may  be  resolved  without 
the  need  to  resort  to  expensive  wind  tunnel  testing.  It  is  therefore  imperative  that  the 
Computational  Fluid  Dynamics  (CFD)  researcher  provide  the  means  to  obtain  accurate, 
timely  solutions  to  these  types  of  flow  fields  if  the  development  of  the  next  generation  of 
aerospace  vehicles  is  to  take  place.  To  ensure  their  long-term  utility,  the  methodologies 
thus  developed  must  be  flexible  enough  to  accurately  rraolve  the  flowfields  over  a  large 
variety  of  geometries  and  under  varying  atmospheric  conditions. 

To  increase  the  confidence  level  in  computational  output  for  arbitrary  geometries 
under  varying  atmospheric  conditions,  steps  must  be  taken  to  validate  the  driver  code’s 
output  with  existing  experimental  data  on  similar  geometries  subjected  to  similar  condi¬ 
tions.  Thus,  valid  data  sets  from  supersonic  turbulent  flow  surveys  play  a  vital  role  in  the 
development  of  turbulence  models  for  CFD.  Recent  experimental  investigations  by  Miller 
(16),  Dotter  (5),  Luker  (15),  Hale  (9)  and  Huffman  (12)  at  the  Air  Force  Institute  of  Tech- 
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nology  (AFIT)  provide  excellent  test  cases  for  the  development  of  turbulence  models  for 
CFD.  The  aforementioned  tests  were  designed  to  meet  the  rigorous  requirements  of  Set¬ 
tles  and  Dodson  (20)  for  viable  high-speed  turbulence  model  validation  cases,  and  as  such 
should  prove  to  be  invaluable  in  the  development  of  turbulence  closure  models.  Since  the 
stated  goal  of  this  data  collection  is  to  serve  as  a  validation  database  for  computational 
turbulence  research,  generating  a  code  with  which  to  perform  this  research  is  a  logical 
companion  effort. 

The  focus  of  this  eflfort  is  threefold:  first,  add  a  two-equation  turbulence  model  to  an 
existing,  validated,  laminar  CFD  code;  second,  validate  the  code  thus  modified  for  use  in 
turbulence  research  by  comparing  its  output  with  the  experimental  data  mentioned  above; 
and  third,  utilize  this  code  to  evaluate  the  effects  of  pressure  gradient  on  compressible, 
turbulent,  boundary  layer  flowfields. 

l.S  The  Closure  Problem 

Turbulence  is  a  three-dimensional  time-dependent  phenomenon,  manifesting  itself  in 
fluid  flows  in  which  the  ratio  of  inertial  forces  to  viscous  forces  is  sufiiciently  large.  In 
1937,  Taylor  and  von  Karman  (see  Wilcox  (27))  proposed  that  “turbulence  is  an  irregular 
motion  which  in  general  makes  its  appearance  in  fluids,  gaseous  or  liquid,  when  they  flow 
past  solid  surfaces  or  even  when  neighboring  streams  of  the  same  fluid  flow  past  or  over 
one  another.”  In  turbulent  flow,  the  viscous  forces  are  insufficient  to  damp  out  small 
disturbances  in  the  mean  flow,  which  instead  are  amplified,  leading  to  irregular,  turbulent 
motion.  The  structure  of  this  turbulent  motion  may  be  thought  of  as  a  cascade  of  eddies, 
in  which  large  eddies  swirling  through  the  flow  entrain  smaller  eddies,  who  in  turn  entrain 
yet  smaller  eddies.  The  energy  associated  with  these  eddies,  the  turbulent  kinetic  energy, 
is  dissipated  away  in  the  smallest  eddies.  This  cascade  of  eddies  leads  to  greater  diffusivity 
in  turbulent  flow,  resulting  in  increased  skin  friction  and  heat  transfer  rates.  The  largest 
of  these  eddies  are  typically  of  the  same  order  of  magnitude  as  the  length  scale  of  the 
mean  fluid  flow;  the  smallest  are  microscopic  in  size.  Thus,  the  range  of  applicable  length 
scales  in  turbulent  flow  is  extremely  large.  Fortunately,  as  noted  in  Wilcox  (27),  “the  time- 
dependent,  three-dimensional  Navier  Stokes  equation  contains  all  of  the  physics  of  a  given 
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turbulent  flow.”  Notice  that  turbulent  fluid  motion  is  not  excluded;  it  too  must  obey  the 
Navier-Stokes  equations.  That  said,  modeling  a  turbulent  flow  should  be  fundamentally  no 
different  than  modeling  a  laminar  flow;  in  fact  turbulence  will  develop  naturally,  provided 
that  all  the  pertinent  length  scales  are  captured. 

But  just  how  widely  do  these  length  scales  vary?  As  previously  noted,  the  largest 
is  on  the  order  of  the  characteristic  length  of  the  flowfield  being  modeled.  The  length 
scale  associated  with  the  smallest  eddies  is  on  the  order  of  the  Kolmogorov  length  scale, 
7/*.  For  wall  bounded  flows,  the  characteristic  length  of  the  flow  may  be  taken  to  be  the 
boundary  layer  thickness,  6,  The  ratio  of  these  length  scales  has  been  found  to  be  inversely 
proportional  to  the  Reynolds  number  based  on  6  raised  to  the  three-fourths  power,  or 


(1.1) 


Using  this  approximation,  i?*  for  the  flowfields  to  be  studied  in  this  investigation  is  esti¬ 
mated  to  be  roughly  1.22xl0~®.  Assuming  that  flowfield  features  of  this  length  scale  may 
be  accurately  resolved  if  the  grid  spacing  is  equal  to  77*,  computing  the  solution  over  1 
cubic  centimeter  of  the  flowfield  would  require  a  grid  with  5.49x10“  nodes.  With  this  kind 
of  required  node  density,  problems  of  a  realistic  size  quickly  become  unmanageable  with 
present  day  computational  resources.  This  is  the  core  problem  in  turbulence.  The  inability 
to  resolve  all  of  the  pertinent  length  scales,  all  the  effects  of  the  turbulent  eddies,  requires 
that  the  effect  of  these  eddies  be  modeled. 


In  present-day  turbulence  modeling,  a  statistical  approach  is  used  to  solve  the  gov¬ 
erning  equations.  In  both  Reynolds  and  Favre  averaging,  each  flow  variable  is  replaced  by 
the  sum  of  a  time-averaged  and  a  fluctuating  term.  Expanding  and  then  contracting  these 
terms  (See  Appendix  E)  produces  similar  sets  of  differential  equations.  These  forms  of  the 
governing  equations,  the  Reynolds-  and  Favre-averaged  Navier-Stokes  equations,  contain 
additional  fluctuation  terms  required  to  model  the  actions  of  the  turbulent  eddies  and  are 
the  forms  most  frequently  used  in  turbulent  flow  modeling. 
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1.3  “Excuse  me,  your  Boussinesq  is  showing” 

At  the  core  of  eddy  viscosity  turbulence  models  is  the  Boussinesq  eddy  viscosity 
approximation,  which  (27)  “assumes  that  the  principal  axes  of  the  Reynolds  stress  tensor, 
Tij  are  coincident  with  those  of  the  mean  strain- rate  tensor,  5^-,  at  all  points  in  a  turbulent 
flow  ...  [where]  ...  the  constant  of  proportionality  between  and  S^-  is  the  eddy  viscosity, 
fir-”  The  ultimate  goal  of  an  eddy  viscosity  turbulence  closure  model  is  then  to  compute 
fpr  in  such  a  manner  as  to  replicate  accurately  the  effects  of  the  turbulent  eddies  on 
the  flowfield.  Dimensional  arguments  have  been  used  to  show  that  appropriate  length 
(or  time)  and  velocity  scales  are  needed  to  specify  fpp.  Two-equation  turbulence  models 
propose  differential  equations  that  govern  these  length  scales  and  as  such  are  complete 
models;  both  required  scales  are  specified.  The  solution  of  these  differential  equations  in 
terms  of  the  scales  modeled  is  then  used  to  specify  nr- 

The  wide  range  of  velocity  and  length  scales  present  in  a  turbulent  flowfield  neces¬ 
sitates  the  consideration  of  non-local  flowfield  influences  on  the  local  turbulent  quantities. 
Turbulent  kinetic  energy  models  have  been  found  to  be  successful  at  incorporating  these 
effects.  Prandtl  chose  the  kinetic  energy  per  unit  mass  of  the  turbulent  fluctuations,  k,  as 
the  basis  for  his  velocity  scale: 


As  described  in  Appendix  E,  the  primed  velocities  in  Equation  1.2  represent  the  Reynolds 
fluctuating  portion  of  the  instantaneous  velocities  and  overbars  represent  time  averaging. 
The  link  between  k  and  nr  lies  in  the  trace  of  the  Reynolds  stress  tensor  (see  Appendix  E), 
which  may  be  expressed  in  terms  of  Favre  fluctuating  components  as 

r7i  =  =  -‘ipk  (1.3) 

The  full  form  of  the  Reynolds  stress  tensor  is  then  given  by 

rj  =  2».  (s.;  -  -  lpkS„  (1.4) 


1-4 


! 

I 


where  the  second  term  ensures  that  the  trace  of  the  Reynolds  stress  tensor  is  correct  for 
incompressible  flows,  where  5,-,  =0.  In  a  similar  fashion,  the  differential  equation  for 
k  may  be  obtained  by  taking  the  trace  of  the  Reynolds  stress  equation  (27).  With  the 
velocity  scale  thus  specified,  a  length  scale  is  needed. 

In  1942,  Kolmogorov  (see  Wilcox  (27))  proposed  a  two-equation  turbulence  model 
that  used  the  specific  dissipation  of  turbulent  kinetic  energy,  u,  as  its  second  scale.  Al¬ 
though  current  trends  lean  towards  viewing  u  as  the  dissipation  of  k  per  unit  k,  oj  itself 
is  a  non-physicaJ  quantity  and  as  such,  its  exact  definition  is  open  for  discussion.  In  justi¬ 
fying  the  selection  of  a;  as  the  second  scale,  Wilcox  (27),  a  strong  proponent  of  the  k  —  u 
turbulence  model,  argues  that  “While  the  actual  process  of  dissipation  takes  place  in  the 
smallest  eddies,  the  rate  of  dissipation  is  the  rate  of  transfer  of  turbulence  kinetic  energy 
to  the  smallest  eddies.  Hence,  this  rate  is  set  by  the  properties  of  the  large  eddies,  and 
thus  scales  with  k  and  /  . . .  ”  Thus,  the  two  scales  modeled  in  the  k  —  u  turbulence  model 
are  consistent;  both  concern  themselves  with  the  actions  of  the  large  eddies  (not  to  be 
confused  with  Large  Eddy  Simulation,  another  turbulence  modeling  method,  which  will 
not  be  discussed  here).  Physical  arguments  have  been  utilized  to  determine  the  form  of  the 
governing  equation  for  w,  which  closely  resembles  the  governing  equation  for  k  in  that  con¬ 
vection,  diffusion,  production,  and  dissipation  terms  are  all  present.  Though  the  arguments 
used  to  determine  the  governing  equation  for  u  are  certainly  not  rigorous  mathematically, 
in  turbulence  modeling  the  most  successful  models  are  frequently  developed  with  more 
emphasis  on  dimensional  analysis  and  physical  reasoning  than  mathematical  rigor. 

The  more  popular  k  —  t  turbulence  model  attempts  to  directly  model  the  dissipation, 
€,  occurring  in  the  smallest  eddies.  The  governing  equation  for  e  is  obtained  by  taking  a 
moment  of  the  Navier-Stokes  equations  and  involves  several  additional  double  and  triple 
correlations.  Since  this  equation  is  obtained  with  substantially  more  rigor  than  the  equation 
for  w,  the  bulk  of  two-equation  turbulent  flow  modeling  has  been  done  with  this  model. 
The  model  is  not  without  flaws,  however.  Closing  the  equation  for  c  requires  estimations 
of  the  aforementioned  correlations,  which  traditionally  have  been  very  difficult  to  measure 
experimentally.  Additionally,  this  difficulty  aside,  there  is  concern  that  even  if  it  is  modeled 
exactly,  c  might  not  be  the  appropriate  length  scale,  since  it  is  the  scale  associated  with 
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the  dissipating  eddies  and  not  that  of  the  Reynolds-stress  bearing  larger  eddies.  Recall 
that  k  is  the  velocity  scale  corresponding  to  the  large-scale  turbulent  fluctuations. 

1.4  Flow  field  Configuration 

The  computations  performed  in  this  investigation  simulate  the  flow  in  the  AFIT 
supersonic  pressure-vacuum  wind  tunnel,  which  operates  at  a  nominal  freestream  Mach 
number  of  2.9.  In  this  facility,  favorable  and  adverse  pressure  gradients  are  imposed  on  su¬ 
personic  flow  through  the  introduction  of  expansion  and  compression  ramps,  respectively, 
into  the  tunnel  ceiling.  In  four  related  studies  at  AFIT,  experimental  data  have  been 
collected  in  these  pressure  gradient  regions.  The  present  numerical  effort  simulates  these 
flowfields  and  compares  the  data  thus  generated  to  that  obtained  through  the  aforemen¬ 
tioned  experimental  investigations. 

1.5  Summary 

The  generality,  completeness,  and  overall  accuracy  of  two-equation  turbulence  models 
makes  them  attractive.  The  Boussinesq  approximation  makes  their  implementation  into 
existing  laminar  codes  easy.  Due  to  questions  regarding  the  basic  formulation  of  the 
k  —  e  model  and  reports  of  exceptional  accuracy  obtained  with  the  k  —  u;  model  (see 
Wilcox  (26)  (27),  Liu  and  Zheng  (14),  Coakley  and  Huang  (4),  and  Morrison  (18)),  the 
focus  of  this  computational  effort  is  on  the  k  —  u  model.  Several  different  versions  of  the 
k  —  Lj  model  exist;  certainly  Wilcox’s  (27)  has  been  tested  more  than  any  other.  For  this 
investigation,  however,  the  model  presented  by  Coakley  and  Huang  (4)  is  used.  The  reason 
behind  this  selection  is  that  the  model  itself  is  fundamentally  the  same  as  Wilcox’s,  but 
the  algorithm  development  presented  in  (4)  is  such  that  only  minor  changes  are  required  to 
move  from  the  k  —  u  model  to  a  A:  —  e  model  or  almost  any  other  two-equation  turbulence 
model.  Thus,  the  final  source  code  will  serve  as  a  valuable  research  tool  for  turbulence 
model  development,  facilitating  direct  comparison  between  results  obtained  with  different 
turbulence  models  using  the  same  basic  flow  solver. 

The  remainder  of  this  thesis  is  structured  as  follows:  A  description  of  the  problem 
geometries  is  presented  in  Chapter  2.  Details  of  the  governing  equations  and  development 
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of  the  numeric  scheme  used  in  the  mean  flow  and  turbulence  modeling  portions  of  the 
source  code  are  provided  in  Chapter  3.  Chapter  four  presents  results  of  a  grid  refinement 
study  conducted  to  ensure  grid  independence  in  the  solutions.  Results  of  the  present  study 
are  compared  and  contrasted  with  experimental  data  collected  by  Miller  (16),  Dotter  (5) 
Luker  (15),  and  Hale  (9)  in  Chapter  5.  Conclusions  and  recommendations,  including 
proposed  coding  changes  are  given  in  Chapter  6. 
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//.  Problem  Descriptions 


2.1  AFIT  Mach  3  Blowdown  Wind  Tunnel 

The  computations  presented  in  this  effort  simulate  the  flowfield  within  the  AFIT 
Supersonic  Pressure- Vacuum  Wind  Tunnel.  This  wind  tunnel  was  used  in  previous  inves¬ 
tigations  by  Dotter  (5)  and  Miller  (16),  who  experimentally  examined  Mach  2.9  boundary 
layer  flows  using  hot-film  anemometry.  Luker  (15)  and  Hale  (9)  considered  the  same 
flowfields  in  their  investigation,  employing  Laser  Doppler  Velocimetry  for  flowfield  charac¬ 
terization.  In  the  latter  set  of  experiments,  the  settling  chamber  pressure  and  temperature 
were  2.1263x10®  Pa  and  294  K,  respectively,  producing  a  freestream  Mach  number  at  the 
nozzle  exit  of  nominally  2.85  and  a  freestream  Reynolds  number  Rg/m  =  1.75x10^.  All 
flow  stations  referenced  in  this  investigation  are  specified  in  centimeters  using  a  Cartesian 
coordinate  system  with  the  origin  set  at  the  nozzle  throat.  As  such.  Flow  Station  (FS)  71.5 
corresponds  to  a  location  71.5  cm  aft  of  the  nozzle  throat.  For  convenience,  computational 
meshes  used  in  this  investigation  are  inverted;  the  wind  tunnel  ceiling,  which  contains  the 
compression  (Adverse  Pressure  Gradient)  and  expansion  (Favorable  Pressure  Gradient) 
ramps  in  the  physical  domain,  is  the  floor  of  the  computational  domain. 

2.2  Nozzle 

Figure  2.1  is  a  cross  section  of  the  nozzle  section  used  in  the  AFIT  Supersonic 
Pressure- Vacuum  Wind  Tunnel.  This  contour  was  designed  by  Dr  Rodney  Bowersox  using 
the  method  of  characteristics  with  a  boundary  condition  correction.  The  coordinates  of 
the  contour  of  this  nozzle  may  be  found  in  Appendix  F. 

2.3  Zero  Pressure  Gradient  (ZPG)  Test  Sections 

Experimental  data  has  been  collected  in  this  test  section,  effectively  a  zero-pressure 
gradient  (ZPG)  channel  flow,  at  FS  44.0  by  Miller  (16)  and  Luker  (15).  Additional  ZPG 
data  was  collected  by  Luker  (15)  in  a  similar  channel  flow  configuration  at  FS  71.5  with 
the  pressure  gradient  test  sections  removed. 


AFIT  Supersonic  Pressure-Vacuum  Wind  Tunnel  Nozzle 


Figure  2.1  Nozzle  Cross  Section 


2.4  Favorable  Pressure  Gradient  (FPG)  Test  Section 

The  tunnel  contour  in  this  test  section  is  presented  in  Figure  2.2  and  is  given  by  the 
third-order  polynomial 


Y  =  -0.2078  4-  0.0897X  -  0.009476a:^  -I-  0.00003598a;®  (2.1) 


Figure  2.2  Favorable  Pressure  Gradient  Test  Section  Cross  Section 

where  x  is  measured  in  cm  from  the  beginning  of  the  test  section  (60.0  cm).  Experimental 
data  has  been  acquired  by  Miller  (16)  and  Luker  (15)  at  FS  71.5  in  a  region  characterized 
by  a  mild  favorable  pressure  gradient.  This  data  plane  is  presented  as  the  dashed  line  in 
Figure  2.2  and  is  locally  normal  to  the  lower  surface.  All  flowfield  quantiti^  at  this  plane 
are  represented  in  terms  of  a  body  normal  coordinate  system,  defined  by  the  dashed  line 
and  the  tangent  line  at  FS  71.5. 
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2.5  Adverse  Pressure  Gradient  (APG)  Test  Section 

Figure  2.4  presents  the  tunnel  contour  in  this  test  section,  which  is  given  by  the 
third-order  polynomial 


Y  =  1.186  -  0.5410a;  +  0.07478x^  -  0.0028x^  (2.2) 

where  x  is  once  again  measured  in  cm  from  the  beginning  of  the  test  section  (60.0  cm). 
Experimental  data  has  been  collected  in  this  test  section  at  FS  71.0  by  both  Dotter  (5) 
and  Hale  (9).  During  the  course  of  this  investigation,  however,  it  became  obvious  that  FS 
71.0  lies  in  a  region  whose  pressure  gradient  is  even  more  favorable  than  that  referred  to  in 
Section  2.4,  as  may  be  seen  in  Figure  2.3.  This  figure  shows  the  computed  static  pressure 
along  the  surface  of  the  contours  of  both  test  sections;  indicating  that  the  pressure  gradient 
at  FS  71.0  of  the  APG  test  section  is  indeed  favorable.  Using  this  figure  to  guide  further 
experiments,  Hale  (9)  obtained  additional  data  at  FS  68.0,  well  within  the  APG  region. 


Figure  2.3  Surface  Pressure  in  Test  Sections 

The  planes  of  data  examined  experimentally  are  shown  by  the  dotted  and  dashed  lines  in 
Figure  2.4.  Once  again,  a  body  normal  coordinate  system  is  used  for  data  output. 
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Figure  2.4  Adverse  Pressure  Gradient  Test  Section  Cross  Section 
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III.  Numerical  Formulation 


3.1  Governing  Equations 

Within  the  field  of  fluid  mechanics,  conservation  of  mass,  momentum,  and  energy 
may  be  expressed  through  the  Navier-Stokes  equations,  which  may  be  written  (in  Cartesian 
coordinates  for  two-dimensional  flow)  as: 


Ut  Fg  Gy  =  0 


(3.1) 


The  solution  vector  is  given  by 


( " '' 
v= 

pv 

and  the  flux  vectors  in  the  x  and  y  coordinate  directions  are 


F  = 


pu'^  +  p- 

pUV  -  Txy 

\  {Et  -I-  p)u  -  «r„  -  vT^y  -I-  g,  / 

^  pv 

pUV  -  T^y 

PV^+P-  Tyy 

\  {Et  -1-  P)t7  -  UT^y  -  VTyy  ^Qy  ) 

In  the  above  formulation,  Et  is  the  total  energy  per  unit  volume  and  is  given  by 


G  = 


Et  =  p 


(3.2) 


(3.3) 


(3.4) 


(3.5) 


where  e  is  the  internal  energy  per  unit  mass.  The  viscous  stresses  and  heat  fluxes  are  given 

by 


2  dut' 


(3.6) 


and 


q.  =  (-«T),  9,  =  i-KT\  (3.7) 

Equation  3.1  is  a  system  of  four  equations  in  eight  unknowns,  those  being  p,  u,  v,  p,  e, 
/i,  T,  and  k.  Thus,  four  additional  relations  are  required  to  close  the  system.  Three  of 
these  come  from  Sutherland’s  law  and  the  jissumptions  of  thermally  and  calorically  perfect 
gases: 


p  =  pRT 


(3.8) 


e  =  c„T 


(3.9) 


rp3/2 

r  +  C2 


(3.10) 


where  R  is  the  gas  constant,  c„  is  the  specific  heat  at  constant  volume,  and  Ci  and  are 
constant  for  a  given  gas.  Assuming  that  the  specific  heat  at  constant  pressure,  Cj,,  and  the 
Prandtl  number,  Pr,  are  constant  for  the  flows  being  considered,  the  final  relation  comes 
from  the  definition  of  the  Prandtl  number: 


K 


(3.11) 
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3.2  Algorithm  Development 

This  section  describes  the  algorithm  used  within  Flat  Plate  Finite  Volume  (FPFV), 
the  driver  code  utilized  and  modified  in  this  investigation.  FPFV  operates  on  the  two- 
dimensional  Navier-Stokes  equations  in  general  curvilinear  coordinates  using  a  cell-centered, 
finite  volume  approach.  This  code  was  originally  developed  and  validated  for  high-speed 
viscous  flows  by  Dr.  Datta  Gaitonde  (8)  (6). 

The  form  of  the  governing  equations  given  in  the  preceding  section  is  not  amenable 
to  solution  on  uneven  computational  meshes.  Transforming  the  uneven  physical  domain 
{{x,y)  space)  into  an  even  computational  domain  ((^,»?)  space)  allows  the  equations  to 
be  more  readily  solved.  As  the  physical  space  is  transformed  into  computational  space, 
so  must  the  governing  equations  be  transformed.  The  set  of  equations  thus  obtained  is 
referred  to  as  the  strong  conservation  form  of  the  governing  equations,  the  derivation  of 
which  may  be  found  in  Appendix  A. 

Beginning  with  the  strong  conservation  form  of  the  governing  equations 


=0 


(3.12) 


the  dissipative  viscous  terms  are  separated  from  the  convective  terms  to  yield  the  form 


where 


1  dU  dFi„„  dGin,  dF,  dG, 
J  dt  '''  '''  dT}  ^  dr] 


(3.13) 


-1-  pV^y 

(PW*  -I-  P)^,  -f-  pUV^y 

pnv^x  +  (pv^  -1-  p)|y 

{Et  +  p)ti6  +  {Et  +  p)v^y 


(3.14) 
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0 


F,= 


7 


~Txy^x  ~  '^yy^y 

{-UTxx  -  VT^y  +  +  (-«T*y  -  VTyy  +  9y)4 


pUTj^  +  pVT)y 

{pv}  +  p)r}^  +  puvpy 
PUVT]^  +  {pv^  +  p)riy 
(Et  +  p)up^  +  (Et  +  p)vT}y 


(3.15) 


(3.16) 


0 

~TxxPx  '^xy'Hy 
‘^xyVx  "^yyVy 

{-UTxx  -  VT^y  +  9,)t/,  +  (-1iT,y  -  VTyy  +  ?„)»?, 


(3.17) 


Using  A  to  represent  spatial  finite  differences,  Equation  3.13  may  be  further  re-written  as 


where 


SU  AF,„,  AGin.  AF,  AG,  _ 
J6t^  A^  At/  A^  At/ 


(3.18) 


SU  =  -  U”  (3.19) 

where  f7”  represents  the  solution  at  the  current  time  step  and  represents  the  solution 
at  the  next  time  step.  Since  the  computational  domain  is  uniform  with  A^  =  At/  =  1.0, 
Equation  3.18  may  be  re-written  yet  again  as 


+  ^KiJ  +  = 


(3.20) 


Solving  for  the  update  to  the  solution  SU  based  on  the  current  state  of  the  flux  vectors 
yields 
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(3.21) 


SU,j  =  V^f' 


[a7?„,,,,  +  AG?. , , + ap; 


In  the  present  computational  effort,  this  equation  is  solved  using  a  cell- centered,  finite- 
volume  approach  in  which  the  conserved  quantities  (p,  pu,  pv,  Et)  are  stored  at  cell  centers. 
Note  that  the  formulation  in  Equation  3.21  may  be  obtained  through  the  use  of  a  classical 
finite-volume  approach  (10)  on  a  uniform  computational  mesh,  as  can  be  found  in  Appendix 
B. 


RolPointt 


SoUBoMdiiy 


OhoAMolf 


j-2 


Figure  3.1  Cell  Configuration  Along  a  Solid  Wall  Boundary 


3.2.1  Boundary  Conditions.  As  the  name  suggests,  cell-centered,  finite- volume 
codes  store  no  flowfield  information  at  the  cell  interfaces.  Specification  of  boundary  condi¬ 
tions  poses  special  difficulties  for  these  codes  since  cell  interfaces,  not  cell  centers,  lie  on  the 
boundaries  of  the  physical  and  computational  domains.  To  address  this  difficulty,  “ghost 
points”  (see  Figure  3.1)  may  be  created  outside  the  physical  domain  and  their  values  set 
to  enforce  the  condition  at  the  interface  (surface).  In  this  effort,  the  physical  domain  is 
broken  into  a  mesh  consisting  of  (il  -  2)  x  (jl  -  2)  cells  with  nodes  2  through  (il  -  1)  in 
the  ^  direction  (2  through  (jl  -  1)  in  the  rj  direction)  treated  as  interior  nodes.  The  ghost 
points  are  assigned  indices  of  1  and  il  in  the  ^  direction  and  1  and  jl  in  the  rj  direction. 

The  boundary  conditions  imposed  at  a  non-moving  impermeable  wall  may  be  ex¬ 
pressed  in  the  following  manner:  No  slip  and  no  penetration  are  specified  by  setting  u 
and  V  at  the  wall  equal  to  zero.  This  is  accomplished  by  specifying  that  the  ghost  point 
velocities  are  equal  in  magnitude  and  opposite  in  direction  to  the  corresponding  velocities 
inside  the  physical  domain: 
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u  UaH  =  0  :  Ui,i  =  -Ui,2 

V  \wall  =  0  :  Vi,l  =  -V,-,2 

In  specifying  an  isothermal  wall  boundary  condition,  the  internal  energy  at  the  ghost  point 
is  extrapolated  using  the  specified  wall  temperature  and  the  internal  energy  at  the  first 
node  inside  the  physical  domain.  The  isobaric  wall  boundary  condition,  ||  =  0,  requires 
that  the  static  pressure  at  the  ghost  point  be  equal  to  the  static  pressure  just  inside  the 
physical  domain: 


T  luioH  —  'I'wall  •  ®f,l  —  2c^T'u,a/|  e,-_2 

UaH  =  0  :  Pi  I  =  Pi  2 

Other  pertinent  flowfield  properties  {p,pu,pv,Et,p)  are  obtained  from  these  specified  val¬ 
ues.  Inflow  boundary  conditions  are  specified  by  explicitly  setting  flowfield  properties  at 
the  upstream  ghost  points: 


Plj  = 

Poo 

«lj  = 

t*0O 

Vlj  = 

Woo 

elj  = 

CvToo 

Plj  = 

Poo 

where  once  s^aJn,  pu,  pv,  Et,  and  p  are  specified  through  the  use  of  these  values  with  the 
governing  equations  listed  in  Section  3.1.  Details  concerning  the  specific  implementation  of 
the  upstream  boundary  conditions  are  presented  in  Section  4.6.  Downstream  ghost  point 
values  are  extrapolated  using  the  values  at  the  last  two  interior  nodes  to  define  a  rate  of 
change  of  each  flow  property,  which  is  then  used  to  set  flowfield  quantities  at  the  ghost 
points: 
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PilJ  —  U,'j_2j- 

pUiij  =  2pUi/_ij  - 

Wj.i  =  'i'P'Oil-lj  -  pVil-2J 
Pu,j  —  ^Pit-lJ  ~  Pit-2J 


^il,i  = 
VilJ  = 

^il,j  = 

Et,ilJ  = 


3.2.2  Implementation.  The  flowfield  is  initialized  using  either  a  uniform  or  one¬ 
dimensional  approximation  of  initial  flowfield  properties.  The  next  step,  the  first  in  an 
iterative  series  of  steps,  is  to  scan  the  flowfield  to  establish  either  a  local  or  global  St  based 
on  a  specified  maximum  CFL  number,  the  magnitude  of  the  mean  flow  eigenvalues,  and 
the  cell  size.  After  St  (potentially  Stij)  has  been  determined,  computational  sweeps  on  the 
interior  nodes  begin  by  fixing  j  and  sweeping  2  <i  <{il  —  1).  In  this  sweep,  denoted  the  ^ 
sweep,  the  left  and  right  states  at  the  cell’s  Eastern  interface  are  extrapolated  using  both 
upstream  and  downstream  flow  information.  MINMOD  limiting  (28)  is  used  in  conjunction 
with  flowfield  information  from  the  t,  (i  -f  1),  and  {*  —  1)  cells  to  establish  the  left  state; 
the  right  state  uses  the  same  limiting  logic  with  flowfield  information  from  the  i,  (t  +  1) 
and  (i  2)  cells.  In  this  application,  p,  pu,  pv  and  p  are  extrapolated;  other  flowfield 
properties  are  obtcuned  from  these  extrapolated  quantities.  The  choice  of  the  extrapolated 
vsuriables  has  a  large  impact  on  the  behavior  of  the  algorithm;  extrapolating  these  variables 
provides  the  best  overall  flowfield  behavior.  See  Appendix  C  for  an  explanation  of  the 
extrapolation  and  limiting  procedure.  The  left  and  right  states  are  combined  using  Roe’s 
flux-difference  splitting  (18)  to  obtmn  the  inviscid  flux  at  the  East  interface  of  each  cell. 
Central  differences  about  the  cell  interface  are  used  to  obtain  the  viscous  fluxes  (including 
the  turbulent  contributions  from  pr  and  Kt)  at  the  East  interface,  which  are  added  to  the 
previously  calculated  inviscid  fluxes.  At  this  point  the  contribution  of  the  East  and  West 
fluxes  to  SU  is  calculated  as  shown  below. 
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inv,i-l,j 


-F. 


(3.22) 


The  T)  sweep  begins  by  fixing  i  and  sweeping  2  <  j  <  (jl  —  1).  The  manipulations 
performed  in  the  t)  sweep  are  effectively  identical  to  those  performed  in  the  ^  sweep  ex¬ 
cept  that  the  states  are  established  at  the  North  interface  of  each  cell,  and  flux  boundary 
conditions  are  enforced  at  solid  wall  boundaries.  These  conditions  supplement  the  bound¬ 
ary  conditions  enforced  explicitly  at  the  end  of  each  step.  The  flux  boundary  conditions 
enforced  at  these  solid  wall  boundaries  are 


0 

-PVi 

0 


(3.23) 


reflecting  «  =  v  =  0  at  the  stationary,  impermeable  walls.  At  the  conclusion  of  these 
sweeps,  the  contribution  of  the  North  and  South  fluxes  to  SU  is  calculated  and  added  to 
that  previously  established: 


—  SUij  -f  JijStij 


•inv.ij 


^inv,ij+l  “b  ^v,ij  ^v,ij+l 


] 


Using  the  SU  vector  calculated  above,  the  flowfield  is  updated  by  setting 


(3.24) 


(3.25) 


At  this  point,  the  residual  is  evaluated  to  determine  the  extent  to  which  the  flow  has 
converged.  The  flowfield  residual  is  defined  as 


RESID  = 


J 


^interior  nodes 


iil-2)ijl~2) 


(3.26) 
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Following  this  calculation,  the  mean  flow  boundary  conditions  are  updated  explicitly  as 
described  in  Section  3.2.1. 

If  desired,  the  above  procedure  may  be  executed  twice  at  each  time  step,  allowing 
the  code  to  operate  in  a  predictor-corrector  mode.  After  the  second  set  of  sweeps,  the 
solution  is  updated  by  setting 


Vt. 


(3.27) 


with  the  residual  evaluated  in  exactly  the  same  manner  as  before.  If  only  one  sweep 
through  the  computational  domain  is  desired,  then 


Uij  =  Uij  (3.28) 

If  the  flowfield  under  consideration  is  turbulent,  (i/r  and  kt  are  updated  at  this 
point  using  either  the  Baldwin- Lomax  or  k  -  u  turbulence  model.  Descriptions  of  these 
turbulence  models  are  found  in  Sections  3.3.1  and  3.3.2,  respectively.  When  the  ratio  of 
the  residual  calculated  in  Equation  3.26  to  the  first  residual  is  below  a  specified  threshold, 
calculations  are  halted  and  the  flowfield  is  deemed  converged. 

3.3  Turbulence  Modeling 

Turbulence  modeling  in  FPFV  is  accomplished  through  the  use  of  either  the  Baldwin- 
Lomax  algebraic  turbulence  model  or  the  k  —u  two  equation  turbulence  model.  Both  are 
eddy  viscosity  models,  making  their  incorporation  into  existing  laminar  codes  a  relatively 
straightforward  exercise.  In  eddy  viscosity  turbulence  models,  the  governing  equations  are 
modified  through  the  addition  of  eddy  viscosity  and  eddy  thermal  conductivity  coefficients 
to  the  already  existing  thermophysical  viscosity  and  thermal  conductivity  coefficients,  as 
seen  below 


fi  =  pi  +  fir  k  =  k  +  kt 


(3.29) 


3-9 


where 


KT’  =: 


CpMr 

Ptt 


(3.30) 


and  Ptt  is  the  turbulent  Prandtl  number,  considered  to  be  a  constant  for  a  given  fluid.  As 
is  common  practice  for  air  at  the  temperatures  considered  in  this  study  (see  Coakley  (4), 
Wilcox  (27),  and  White  (25)),  a  value  of  0.90  has  been  assumed. 


3.3.1  Baldwin  -  Lomax  Turbulence  Model.  The  choice  of  an  algebraic  turbulence 
model  as  the  “other”  model  to  be  implemented  in  this  effort  is  based  on  the  well  docu¬ 
mented  features  of  this  flowfield.  The  work  of  Miller  (16)  and  Dotter  (5)  established  that 
the  test  sections  are  separation-free,  meaning  that  simpler  algebraic  turbulence  models 
should  do  a  fair  job  at  computing  the  flowfield.  The  results  from  computations  performed 
with  this  model  may  serve  then  as  an  additional  reference  point  for  the  validation  of  the 
more  complex  k  —  u  model.  The  Baldwin-Lomax  turbulence  model  is  a  two-layer  eddy 
viscosity  model,  incorporating  different  formulations  for  the  eddy  viscosity  for  the  inner 
and  outer  regions  of  the  boundary  layer.  To  that  end,  the  eddy  viscosity  is  defined  as 


y<ym 

y>ym 


(3.31) 


where  is  the  first  point  away  from  the  wall  at  which  fPr^  =  fir,.  This  model  was 
incorporated  in  accordance  with  the  formulation  described  in  Wilcox  (27).  The  sections 
below  describe  the  formulations  in  each  of  these  regions 


3.3.1. 1  Inner  Region.  The  inner  region  uses  a  basic  mixing  length  model 


PTi  =  P^mir  \^Bl\ 


(3.32) 


where  u  is  the  magnitude  of  the  vorticity  vector  and  is  defined  for  two-dimensional  flow  as 
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^BL  — 


/  dv 


dy) 


1/2 


and  the  mixing  length  is  defined  as 


with 


^mti  =  i^blV 


A+  =  26.0 


(3.33) 


(3.34) 


«+  =  ^ 


(3.35) 


_  I'^viaU 

~  V  pivall 


(3.36) 


/ \| 

■ "  Is; + fe)Li 


(3.37) 


3. 3. 1.2  Outer  Region.  In  the  outer  region,  the  Baldwin-Lomax  model  uses 
the  vorticity  in  the  boundary  layer  to  determine  the  length  scale,  giving 


Mr,  =  PoCcpF^akeFKM  ivi 

\  ^Kleh  / 


(3.38) 


where 


Fkub  (y; ^)  = 


1  +  5.51J 


-1 


(3.39) 


■flafce  =  rnin 


Vniax  Ffnax  1  CxukVma 


ui 


dij 


(3.40) 
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(3.41) 


Fmax  =  l‘^flX-|)] 

and  ymax  is  the  y  location  at  which  F^ax  is  evaluated.  Traditionally,  U^ij  is  the  maximum 
value  of  U  for  boundary  layer  flows.  In  this  application,  Udif  is  chosen  to  be  the  centerline 
velocity.  Additional  closure  coefficients,  taken  from  Wilcox  (27),  are  as  follows 

Ckuh  =  0-3 
=  1.6 
C^t  =  1.0 
a  =  0.0168 
Kbl  =  0.41 


Note  that  =  26  is  “rigorously”  valid  only  for  zero  pressure  gradient  flows;  modifications 
may  be  made  to  enhance  the  performance  of  this  method  in  the  presence  of  strong  pressure 
gradients  (25).  In  this  investigation,  the  pressure  gradients  are  mild  and  the  aforementioned 
modification  has  not  been  made. 

3.S.1.3  Implementation.  The  flowfields  under  consideration  are  treated 
as  horizontally  opposed  flat  plate  flows.  Therefore,  at  each  t-plane  in  the  computational 
domain  two  sweeps  aje  used  to  calculate  pr-  The  first  sweep  begins  at  the  top  wall, 
calculating  pn  from  the  wall  to  the  centerline  from  above.  Following  this  calculation,  pr, 
is  calculated  over  the  same  range  of  nodes  and  compared  to  pt-  at  each  node.  The  stored 
values  for  pr  in  the  upper  half  of  the  computational  domadn  are  then  determined  through 
the  rule  given  in  Equation  3.31.  The  above  procedure  is  repeated  on  the  bottom  half  of 
the  computational  domain,  with  the  sweeps  moving  from  the  lower  wall  to  the  centerline 
from  below.  In  this  application,  y  is  taken  to  be  the  normal  distance  from  the  node  in 
question  to  the  node  on  the  wall  at  that  »-plane.  Calculation  of  pr  progresses  in  this 
manner  through  the  computational  domain  one  t-plane  at  a  time  after  each  mean  flow 
update. 
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s.s.s  k  -  „  Turbulence  Model.  As  described  m  Coatley  and  Huang  (4)  the 
governing  differential  equation  for  the  *  -  c  turbulence  model  may  be  written  as 


^KE  ,  dFKB.in,  ,  dGfCE.in.  .  dF, 

"  I  7i  ■  + 


dt 


dx 


dy 


QGke.v  _  p 
dx  dy  ~  “ 


(3.42) 


with 


T 


du  dv 
dx  ^  dy 


(3.49) 
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(3.50) 


with  closure  coefficients  given  by 

(Tk  =<Tu  =  2.0 

fk  =  L  =  f=^  1.0 

Otk  —  fC'il 

Cti  =  C'ta  =  1.0 


C^i  =  0.555 
Ca,2  =  0.833 


Following  the  solution  of  these  differential  equations,  the  eddy  viscosity  may  be  obtained 
through  the  following  relation 


Ht  =  (3.51) 

with  Cft  =  0.09 

3. 3. 2. 1  Algorithm  Development.  The  transformation  of  the  governing  equa¬ 
tions  for  the  k  —  (jj  turbulence  model  into  strong  conservation  form  is  analogous  to  the 
transformation  of  the  governing  mean  flow  equations.  The  resulting  form  of  the  equations 
is 


1  dUfCE  dFi„„  dF,  dG,  P 

j  at  di  dv  dr)  j 

where 


IP  _  {*J’'KE,i«»+^»GKE,in.  _  >IxfKE.ia>+>7«GKE.>»« 

f  inv  —  - ' - j  - ' -  <J«n«  —  - ' - j  - ' - 

Ip  {«f*KE,«i»e+t»GKE.»i»c  7^  VxFk B .wiic+VtG K B .titc 

V  —  J  —  J 


(3.52) 


(3.53) 


As  in  Equation  3.18,  A  is  used  to  represent  spatial  differentiation,  allowing  Equation  3.52 
to  be  written  in  finite-difference  form  as 


(3.54) 


^UkE  ,  ^Finu  ^Gjnv  ,  A(j„  _  P 

J6t  Atj  AC  At?  ~  y 

Defining  SUke  ^  had  been  done  previously  with  in  Equation  3.19,  multiplying  through 
by  J6t,  taking  advantage  of  the  uniform  and  unitary  the  computational  domain,  and 
evaluating  the  flux  vectors  and  source  term,  P,  at  the  current  time  step.  Equation  3.54 
may  be  re-written  as 


^UKE,i.j  -  AP +  AF „  -I-  AG„  (3.55) 


With  all  flow  and  turbulence  properties  known  at  time  level  n,  updating  the  solution 
involves  evaluation  of  the  above  equation  at  each  cell. 

3. 3. 2. 2  Boundary  Conditions.  Boundary  conditions  in  this  application  are 
implemented  in  accordance  with  the  algorithm  description  given  in  Coakley  and  Huang  (4). 
At  an  impermeable  wall,  the  governing  equations  for  k  and  or  must  satisfy 


Jfc  =  0 


(3.56) 


Wi  =  7.2^  (3.57) 

y? 

In  Equation  3.57,  the  subscript  1  refers  to  the  first  calculated  data  point  above  the  wall  and 
y  is  the  normal  distance  from  the  wall  to  that  data  point.  Inflow  or  freestream  conditions 
must  satisfy  the  somewhat  more  ambiguous  conditions 


\  A*  /oo  /oo 


<  1.0 


(3.58) 


Woo  >  10- 


(3.59) 
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where  X  is  a  characteristic  length  of  the  flow  field,  taken  to  be  1  m  in  the  present  investiga¬ 
tion.  For  the  present  calculations,  discussions  with  Dr  Gaitonde  (7)  lead  these  boundary 
conditions  to  the  form 


*00  =  0.001^-^-  (3.60) 

Woo  =  10^  (3.61) 

3. 3. 2.3  Implementation.  The  procedure  by  which  the  governing  equations 
for  the  k  —  u>  turbulence  model  are  marched  to  convergence  is  essentially  identical  to  that 
used  for  the  mean  flow.  As  in  the  mean  flow  equations,  the  conserved  variables  (pk  and 
pu)  are  extrapolated  to  the  cell  faces  using  MINMOD  limiting  and  the  resulting  fluxes 
established  through  the  use  of  Roe  flux-difference  splitting  (18).  The  viscous  terms  are 
once  ^ain  treated  with  central  differences  and,  as  stated  previously,  the  source  terms  are 
evaluated  at  the  present  time  level. 

The  boundary  conditions  are  treated  in  a  similar  manner.  At  the  solid  wall  bound¬ 
aries,  Equation  3.56  is  enforced  by  setting 

ki,i  =  -*.,2  (3.62) 

Since  the  w  solid  wall  boundary  condition  is  not  actually  enforced  at  the  wall,  but  at  the 
point  immediately  above  the  wall.  Equation  3.57  is  implemented  by  setting 

Wi,2  =  7.2^  (3.63) 

y? 

where  yi  is  the  normal  distance  from  the  wall  to  the  cell  center  of  node  (i,2).  Inflow 
boundary  conditions  are  set  using  the  methodology  outlined  in  Section  4.6.  Downstream 
or  outflow  boundary  conditions  are  set  by  extrapolating  k  and  u  using  the  previous  two 
upstream  points  in  a  manner  reminiscent  of  the  mean  flow  calculations.  Subsequent  to  the 
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evaluation  of  the  boundary  conditions  on  k  and  w,  fj^  is  evaluated  through  the  functional 
form  given  in  equation  3.51. 


IV.  Grid  Refinement  and  Flat-Plate  Validation 

A  grid  study  was  conducted  to  ensure  that  the  computed  results  are  not  grid  de¬ 
pendent.  For  this  study,  a  series  of  9  grids  were  examined  to  determine  the  effects  of 
grid  properties  on  the  computed  solutions.  The  properties  of  the  grids  studied  may  be 
found  in  Table  4.1.  The  grids  itemized  in  this  table  begin  at  Flow  Station  (FS)  27.457, 
corresponding  to  the  nozzle  exit  plane,  and  extend  to  FS  73.0,  aft  of  the  favorable  pressure 
graxlient  test  (FPG)  section.  Flowfield  properties  are  examined  both  upstream  at  FS  44.0 
and  in  the  FPG  test  section  at  FS  71.5.  For  the  purposes  of  this  grid  study,  the  upstream 
boundary  conditions  reflect  a  Af  =  2.9,  Rtfm  =  1.5x10^  uniform  flow  with  zero  boundary 
layer  thickness  and  T^aii  =  To  =  294K  impinging  on  the  leading  edge  of  the  grid.  These 
conditions  reflect  the  freestream  conditions  reported  by  Miller  (16)  and  Dotter  (5). 

GRIDGEN2D  revision  8.4.7.1  was  utilized  for  the  grid  generation  task  in  this  effort. 
Node  stretching  normal  to  the  wall  was  accomplished  within  GRJDGEN2D  using  Vinokur’s 
weighting  factor.  Initial  wall  spacings  are  the  same  at  the  North  and  South  grid  boundaries 
and  are  as  given  in  Table  4.1.  Spacing  in  the  streamwise  direction  is  uniform.  Grids  were 
generated  algebraically  using  the  arclength-based  Trans-Finite  Interpolation  (TFI)  option 
and  smoothed  elliptically  using  the  Thomas  and  Middlecoff  background  control  function. 
Transverse  node  spacing  was  held  constant  through  the  smoothing  process,  which  enforced 
orthogonality  at  domain  boundaries  and  utilized  exponential  blending  between  boundary 
and  interior  points.  Grids  were  smoothed  elliptically  until  the  residual  dropped  below 
1.0x10-®. 

The  effects  of  variations  in  streamwise  and  normal  node  count  on  non-dimensional 
velocity  and  turbulent  shear  stress  profiles  and  law  of  the  wall  plots  are  discussed  in  this 
chapter.  An  assessment  of  the  sensitivity  of  the  computed  solution  to  the  initial  wall 
spacing  is  aJso  included. 

4-1  Variation  in  Normal  Node  Count 

The  effect  of  normal  node  count  variation  is  studied  by  doubling  the  normal  node 
count,  utilizing  grids  GSl,  GS2,  and  GS3.  Streamwise  spacing  in  all  three  grids  is  uniform. 
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Table  4.1  Dimensions  -  Grid  Study  Meshes 


Reference 

Size 

^ymid,L 

[m] 

^ymid,R 

[m] 

GSl 

101  X  101 

lESEBI 

GS2 

101 X  151 

t'ni'TTITTTM 

GS3 

101 X  201 

iTTTTMi 

jUAJHUUdl 

1‘ltTi‘ITTt‘il 

GS4 

75  X  151 

iTTTTT'li 

GS5 

151 X  151 

GS6 

201  X  151 

lim^ariBa 

I'nTiTnvM 

GS7 

101  X  151 

GS8 

101 X  151 

fn»{T.^gnigi 

wmmuaM 

GS9 

101  X  151 

^20111 

■iwiKiianfcdl 

esismebi 

and  the  initial  wall  spacing  is  1.5xl0~®m,  calculated  to  provide  values  near  1.0  in  the 
test  sections.  All  figures  referred  to  in  this  section  may  be  found  in  Appendix  G. 

Figures  G.l  and  G.2  show  the  influence  of  normal  node  count  on  the  wall  law  plot, 
which  plots  u'*'  as  a  function  of  y'*’,  where  is  based  on  the  elfective  velocity  as  defined 
by  van  Driest  (25) 


«e//=/  (4.1) 

Jo  \PwaU/ 

Differences  between  the  three  grids  are  very  slight  in  the  inner  region,  which  may  be 
attributed  to  the  fact  that  initial  node  spacing  is  held  constant  for  all  three.  Grid  GSl 
separates  itself  from  the  others  in  the  logarithmic  and  wake  regions,  where  GS2  and  GS3 
continue  to  agree  very  closely.  Nondimensional  velocity  profiles  in  these  two  re^ons  exhibit 
the  same  basic  trends.  Figure  G.3  and  Figure  G.4  are  plots  of  the  nondimensional  velocity 
through  the  boundary  layer.  The  boundary  layer  thickness,  is  the  location  at  which 
u  =  0.99517e,  where  17*  is  the  edge  velocity  at  the  flow  station  of  interest.  Agreement 
between  the  three  solutions  is  good  near  the  wall;  once  again  the  coarser  grid  does  not  do 
well  away  from  the  wall,  where  the  disparity  between  that  solution  and  the  others  reaches 
0.5%.  Comparisons  of  the  nondimensional  turbulent  shear  stress  profiles 
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where  /iy  is  the  eddy  viscosity,  as  defined  in  Chapter  III,  are  shown  in  Figures  G.5  and 
G.6.  Once  again,  differences  between  grids  GS2  and  GS3  are  negligible.  The  figures 
discussed  above  indicate  that  while  101  normal  nodes  is  insufficient  for  accurate  flowfield 
computation,  increasing  the  normsd  node  count  beyond  151  nodes  has  a  negligible  effect  on 
computational  accuracy.  For  this  reason,  the  bulk  of  the  computations  performed  hereafter 
utilize  151  nodes  in  the  normal  direction. 


4.8  Variation  in  Streamwise  Node  Count 

In  this  section  the  effects  of  variations  in  streamwise  node  count  are  studied  by 
examining  computational  output  generated  on  grids  GS2,  GS4,  GS5,  and  GS6.  Uniform 
streamwise  spacing  is  once  again  used  for  the  four  grids  considered,  which  range  in  size 
from  75x151  (Ax  =  .6244cm),  to  201x151  (Ax  =  .2291cm).  Initial  wall  spacing  for  all 
grids  is  once  again  1.5x10"®  m. 

Figures  G.7  and  G.8  are  wall  law  plots  for  the  upstream  flat  plate  and  favorable 
pressure  gradient  test  sections,  respectively.  With  the  exception  of  the  extremely  coarse 
75x151  grid,  axial  node  distribution  seems  to  have  no  effect  on  the  flowfield  properties 
in  the  flat-plate  region.  Figure  G.8,  however,  shows  that  at  FS  71.5  in  the  FPG  test 
section,  axial  node  count  has  a  more  noticeable  effect  on  the  solution  in  the  inner  and 
overlap  regions.  The  difference  disappears  away  from  the  wall  and  thus  seems  to  be  linked 
to  poor  resolution  of  the  curvature  of  the  wall.  The  nondimensional  velocity  profiles  in 
Figure  G.7  and  Figure  G.8  show  little  sensitivity  to  the  streamwise  node  count,  as  do  the 
nondimensional  turbulent  shear  stress  profiles  shown  in  Figures  G.ll  and  Figure  G.12. 

In  order  to  improve  resolution  in  the  curved  wall  region  without  wasting  grid  points 
upstream,  a  non-uniform  node  spacing  scheme  is  adopted  in  the  final  grid  generation  task. 
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4-3  Variation  in  Wall  Spacing 

In  this  section,  computations  performed  on  grids  GS2,  GS7,  GS8,  and  GS9  are  used 
to  determine  the  flowfield  effects  of  initial  wall  spacing  variations.  Uniform  streamwise 
spacing  is  used,  and  the  initial  wall  spacing  is  varied  from  0.5x10"®  m  to  2.0x10"®  m.  The 
effect  of  wall  spacing  on  the  wall  law  plots  may  be  seen  in  Figures  G.13  and  G.14  for  flow 
stations  44.0  and  71.5,  respectively.  These  figures  clearly  show  that  the  initial  wall  spacing 
has  an  effect  on  the  computed  solution.  At  FS  44.0,  the  profiles  generated  using  grid  GS9 
are  appreciably  different  from  the  remainder  of  the  grids.  Moving  to  FS  71.5,  this  difference 
diminishes,  but  does  not  disappear.  Figure  G.15  and  Figure  G.16  show  the  effect  of  wall 
spacing  variation  on  the  nondimensional  velocity  profiles.  Once  again,  differences  at  FS 
44.0  are  greater  than  those  at  FS  71.5,  and  both  curves  show  close  agreement  between  all 
solutions  near  the  wall.  Figures  G.17  and  G.18  serve  only  to  reinforce  the  previous  trends. 

In  considering  initial  wall  spacing  variation  only,  differences  in  the  computed  solutions 
are  greatest  in  the  logarithmic  and  outer  regions.  This  can  be  attributed  to  the  fact  that 
when  the  total  node  count  is  held  constant,  decreasing  the  initial  node  spacing  serves  also 
to  decrease  the  number  of  nodes  that  are  avmlable  for  freestream  resolution.  A  compromise 
must  be  made  between  resolution  of  mean  flow  and  near  wall  effects.  Grids  utilized  in  the 
remainder  of  the  calculations  contained  herein  use  an  initial  wall  spacing  of  1.0x10"®  m, 
providing  for  both  accurate  near-wall  modeling  and  mean  flow  resolution. 

4.4  Final  Grid  Generation 

The  preceding  grid  study  was  used  to  guide  the  final  grid  generation  task  for  this 
project.  All  of  the  grids  begin  with  a  channel  flow  section,  basically  two  parallel  fiat  plates, 
beginning  at  either  FS  0.0  (FPGl,  APGl)  or  FS  44.0  (FPG2,  APG2).  The  streamwise 
node  spacing  in  these  sections  is  uniform  and  given  by  Axj  in  Table  4.2.  Aft  of  this 
section  is  a  4.0  cm  region  where  the  streamwise  spacing  is  compressed  from  Aii  to  Axa 
using  Vinokur’s  weighting  factor.  Axa  is  the  uniform  streamwise  node  spacing  used  in 
the  test  sections,  which  begin  at  FS  64.0  and  extend  to  FS  73.0.  As  was  mentioned  in 
Section  4.2,  node  spacing  in  the  curved  wall  test  sections  is  significantly  tighter  than  that 
used  in  the  upstream  region.  For  grids  FPGl  and  APGl,  the  streamwise  spacing  in  the 
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Table  4.2  Grid  Dimensions  -  Test  Case  Grids 


Reference 

Size 

^  X  q 

^ywall 

[m] 

Axi 

[m] 

FPGl 

219  X  151 

FPG2 

151  X  151 

mQjjsii 

gili*T»Wr>Bg 

APGl 

219  X  151 

IISKSiSM 

APG2 

151  X  151 

FLATl 

151  X  101 

N/A 

test  sections  is  equivalent  to  utilizing  approximately  278  evenly  spaced  nodes  on  the  GS 
series  of  grids.  For  the  FPG2  and  APG2  grids,  this  number  climbs  to  approximately  370 
evenly  spaced  nodes.  This  node  distribution  scheme  enables  increased  resolution  of  the 
curved  wall  sections  without  adversely  impacting  computational  expense.  The  compressed 
region,  including  the  test  sections,  of  Grids  FPGl,  FPG2,  APGl,  and  APG2  are  depicted 
in  Figures  G.19,  G.20,  G.21,  and  G.22,  respectively. 


4.5  Flat  Plate  Validation 

Prior  to  comparisons  with  the  experimental  data  taken  at  AFIT,  validation  runs  using 
the  Baldwin- Lomax  and  k  —  cj  turbulence  models  were  performed  against  Morrison’s  (18) 
flat-plate  test  case  using  the  following  conditions 


M^: 

2.0 

222K 

Ty/all  ‘ 

222K 

Re/m : 

1.0x10 

These  computations  were  performed  on  grid  FLATl,  which  utilized  node  clustering  at  the 
North  and  South  domain  boundaries  and  even  spacing  in  the  streamwise  direction.  This 
nodal  distribution  was  chosen  to  correspond  approximately  to  that  used  in  the  test  section 
calculations.  Results  of  these  runs  may  be  seen  in  Figures  4.1  and  4.2.  Figure  4.1  is  a 
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iH>  vs  y«.  van  Drtest  II  Thaory 


Figure  4.1  y'^  vs  ti+(van  Driest) 


plot  of  vs  I/+,  where  u"*"  is  once  again  based  on  the  effective  velocity  as  defined  by  van 
Driest  (25).  The  correlated  curve  shown  follows  the  rule 


«■•■  =  < 


JL 

.41  “ 


(4.3) 


as  given  in  Morrison  (18).  Also  shown  in  this  plot  is  the  output  from  a  previously  validated 
finite-difference  code  (19)  which  utilizes  Roe  flux  difference  splitting  and  &  k  —  e  two 
equation  turbulence  model.  Agreement  between  the  present  calculations,  the  k  —  e  code, 
and  the  correlated  curve  is  very  close  through  the  inner  region.  In  the  overlap  region, 
both  the  k  —  0}  and  Baldwin-Lomax  models  predict  slightly  lower  velocities  than  the  k  —  e 
calculation.  Note  that  the  correlated  curve  has  not  been  modified  to  account  for  the  overlap 
region  between  the  inner  and  logarithmic  regions.  Both  of  the  two  equation  models  do 
exceptionally  weU  in  the  logarithmic  layer,  where  both  the  magnitude  and  slope  of  the 
correlated  curve  is  replicated  computationally.  In  this  region,  the  Baldwin-Lomax  model 
accurately  predicts  the  slope  of  the  correlated  data,  but  the  magnitude  is  slightly  lower 
than  expected.  Finally,  agreement  between  all  three  models  is  excellent  in  the  wake  region. 


Figure  4.2  Cj  vs  Reg 


Figure  4.2  is  a  plot  of  Cf  vs  ,  where  0  is  the  momentum  thickness  and  is  defined 
(25)  as 


0 


(4.4) 


A  correlated  curve  from  Hopkins  and  Inouye  (11)  is  presented  along  with  the  computational 
results  to  serve  as  a  basis  for  comparison.  All  three  models  accurately  predict  the  skin 
friction  for  iZ*,  >  7500.  At  lower  values  of  Rg,,  the  k  —  €  model  slightly  under-predicts  the 
skin  friction  whereas  the  k  —  tjj  and  Baldwin-Lomax  results  continue  to  agree  exceptionally 
well  with  the  correlated  values. 


4.6  Upstream  Boundary  Conditions 

In  accordance  with  the  recommendations  of  Settles  and  Dodson  (20),  upstream 
boundary  conditions  are  set  using  experimentally  determined  values.  Calculations  per¬ 
formed  in  this  investigation  utilize  two  sets  of  boundary  conditions,  both  based  on  the  FS 
44.0  data  of  Luker  (15).  The  method  by  which  these  conditions  were  set  is  presented  in 
this  section. 
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4-6.1  Freestream  Inflow.  This  boundary  condition  is  used  exclusively  with  grids 
FPGl  and  APGl.  It  assumes  that  an  undisturbed  freestream  =  0)  flow  impinges  on  the 
leading  edge  of  these  grids  at  FS  0.0.  The  freestream  value  of  u  is  set  to  the  edge  velocity 
measured  by  Luker  (15)  at  Flow  Station  (FS)  44.0  and  v  is  assumed  to  be  everywhere  zero. 
The  total  temperature  (To  =  294  K)  and  pressure  (po  =  2.12634x10®  Pa)  are  assumed 
constant,  allowing  direct  computation  of  the  static  temperature  and  pressure  at  this  flow 
station.  The  second  of  these  assumptions  (po  =  constant)  is  used  with  the  realization  that 
it  is  technically  inaccurate  due  to  the  presence  of  weak  shocks  in  the  wind  tunnel  caused 
by  test  section  misalignment. 


=  «oo  —  605.1  m/s 

(4.5) 

To  -  —  =  111.8  K 

2cp 

(4.6) 

.  (0^=7210  Pa 

(4.7) 

The  viscosity  and  density  are  calculated  using  Sutherland’s  Law  3.10  and  the  Equation  of 
State  3.8,  respectively.  The  constants  utilized  in  the  above  formulations  are 

R  =  287.0 
Cp  =  1005.0 
7  =  1.40 
Cl  =  1.458x10-® 

Cj  =  110.4 


J 

EglT 

J 

EpT 

ILs 

m> 

K 


Boundary  conditions  on  k  and  u  are  specified  using  Equations  3.60  and  3.61,  respectively, 
using  the  thermophysical  conditions  calculated  above. 

4-6.2  Experimental  Profile  Inflow.  This  boundary  condition  set  is  used  exclu¬ 
sively  with  grids  FPG2  and  APG2,  approximating  the  flowfield  conditions  measured  by 
Luker  (15)  at  FS  44.0.  The  experimental  velocity  profile  (with  the  additional  specification 
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that  t;  =  0)  is  used  with  the  assumption  of  constant  To  to  approximate  the  temperature 
distribution  through  the  boundary  layer.  The  static  pressure  is  assumed  to  be  constant  and 
equal  to  the  value  calculated  at  the  mid-plane  using  the  methodology  introduced  above. 
From  these  values,  the  density  distribution  through  the  boundary  layer  is  established. 

Boundary  conditions  on  k  are  specified  from  fluctuating  velocity  information  provided 
in  the  experimental  data  set,  where  the  contribution  of  w'  is  assumed  to  be  equal  to  v', 
reducing  Equation  1.2  to 


k  =  ^(u'^  +  2v'^)  (4.8) 

The  velocity  profile  is  differentiated  to  obtain  |^,  which  is  combined  with  the  experi¬ 
mentally  determined  turbulent  shear  stress  profile  to  specify  an  effective  (Pr  distribution. 
Rearranging  Equation  3.51,  k  and  nr  are  used  to  generate  the  u  profile  at  this  station.  The 
k  and  u  profiles  thus  generated  may  be  found  in  Figures  4.3(a)  and  4.3(b),  respectively. 


(a)  (b) 

Figure  4.3  k  and  w  Experimental  Inflow  Profiles 
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V.  Results  and  Discussion 


This  chapter  presents  the  results  of  the  present  computational  effort,  beginning  with 
a  validation  of  the  computed  results  through  a  comparison  with  the  experimental  data  of 
Miller  (16),  Dotter  (5),  Luker  (15)  and  Hale  (9).  Comparisons  are  made  in  the  Zero  Pres¬ 
sure  Gradient  (ZPG),  Favorable  Pressure  Gradient  (FPG),  and  Adverse  Pressure  Gradient 
(APG)  test  sections.  Through  this  analysis,  the  effects  of  pressure  gradient  are  addressed. 
Also  discussed  are  differences  between  the  Baldwin-Lomax  and  k  —  u  solutions,  relative 
computational  expense  for  both  models,  and  the  effect  of  upstream  boundary  condition 
specification  on  the  downstream  solutions. 

5.1  Validation 

This  section  presents  the  results  of  computations  performed  using  the  k  —  u  turbu¬ 
lence  model.  Computations  utilize  grids  FPG2  and  APG2,  and  in  accordance  with  the 
recommendations  of  Settles  and  Dodson  (20),  the  experimental  profile  inflow  conditions 
are  used  unless  otherwise  specified.  Overall  agreement  with  experimental  data  is  found  to 
be  excellent. 


5.1.1  ZPG  Test  Section.  The  inflow  plane  for  grids  FPG2  and  APG2  is  at  Flow 
Station  (FS)  44.0;  as  such,  direct  comparison  to  the  experimental  data  collected  at  FS 
44.0  is  meaningless.  In  the  figures  that  follow,  computed  results  at  FS  60.0,  upstream 
of  the  APG/FPG  test  sections,  is  compared  to  experimental  data  collected  at  FS  44.0 
(Miller  (16)  and  Dotter  (5))  and  FS  71.5  (ZPG  data  from  Luker  (15)  and  Hale  (9)).  This 
comparison  assumes  that  the  boundary  layer  at  all  three  stations  is  fully  developed  and 
therefore  velocity  and  shear  stress  profiles  at  these  stations  should  be  similar  when  plotted 
against  y/6. 

Figure  5.1  is  a  comparison  between  the  present  computed  Mach  number  profile  at 
FS  60.0  and  Miller’s  (16)  experimental  data,  collected  at  FS  44.0.  The  Mach  number  is 
based  on  the  magnitude  of  the  velocity  vector  and  the  local  static  temperature 
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(5.1) 


M  = 


~7W~ 


Agreement  between  the  computational  and  experimental  profiles  is  very  good,  with 
a  slight  (magnitude  is  less  than  1.1%)  mismatch  in  the  freestream  that  may  potentially  be 
related  to  either  a  mismatch  in  experimental  testing  conditions  or  the  physical  separation 
between  the  experimental  and  computational  flow  stations  used  to  generate  the  plot.  With 
regard  to  the  former,  recall  that  the  upstream  boundary  conditions  for  this  calculation  are 
set  using  Luker’s  (15)  FS  44.0  data.  Using  this  inflow  information,  the  freestream  Reynolds 
number  per  meter  is  approximately  1.75x10^,  whereas  Miller  (16)  reported  a  freestream 
Reynolds  number  per  meter  of  1.50x10^  for  his  measurements. 


Figure  5.1  Mach  Number  vs  vISm,  ZPG,  FS  60.0/44.0 

Figure  5.2  is  a  plot  of  the  nondimensional  turbulent  shear  stress  (see  Equation  4.2) 
versus  vISm-  Shown  in  this  plot  are  the  present  computed  results  and  the  “incompressible” 
component  of  the  Reynolds  stress  tensor  (see  Appendix  E)  as  reported  by  Miller  (16)  for 
the  flow  stations  indicated.  Agreement  here  is  quite  good,  with  the  character  of  the  exper¬ 
imental  data  being  accurately  reproduced  computationally.  Vertical  shifting  of  the  data 
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may  be  attributed  to  the  methodology  used  to  determine  the  boundary  layer  thickness,  Sm. 
The  boundary  layer  thickness  utilized  in  this  figure  and  all  those  generated  by  Miller  (16) 
and  Dotter  (5)  is  ^.995^5  or  the  value  of  y  at  which  M  =  0.995Me.  This  process  hinges 
on  the  accurate  determination  of  the  edge  Mach  number,  an  evaluation  made  extremely 
difficult  by  scatter  in  the  experimental  data.  In  retrospect,  defining  the  boundary  layer 
edge  as  ^.99  would  have  significantly  simplified  this  task  by  reducing  the  effect  of  the  data 
scatter.  With  this  in  mind,  computational  agreement  with  the  hot  film  data  of  MiUer  (16) 
at  this  flow  station  is  considered  excellent. 


Figure  5.2  Nondimensional  Turbulent  Shear  Stress  vs  ZPG,  FS  60.0/44.0 

Figure  5.3  plots  the  nondimensionalized  u  component  of  velocity  versus  where 
Su  is  the  value  of  y  at  which  tt  =  0.995uc.  This  figure  and  all  figures  contadning  the 
Laser  Doppler  Velocimetry  (LDV)  data  of  Luker  (15)  and  Hale  (9)  use  this  definition  of 
the  boundary  layer  thickness.  Agreement  between  the  present  calculations  and  the  LDV 
data  is  excellent.  Both  the  magnitude  and  slope  of  the  experimental  data  are  replicated 
computationally.  Similarly  impressive  results  are  presented  in  Figure  5.4,  where  the  nondi¬ 
mensional  turbulent  shear  stress  is  plotted  versus  Once  again,  both  the  shape  and 

magnitude  of  the  experimental  curve  are  captured  in  the  computational  results.  Inspection 
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of  Figure  5.5  reveals  that  the  computed  turbulent  kinetic  energy  (k)  profile  closely  matches 
the  experimental  data  throughout  the  boundary  layer. 


Figure  5.3  u/«e  vs  y/6u,  ZPG,  FS  60.0/71.5 


5.1.2  FPG  Test  Section.  Profiles  presented  in  this  section  are  taken  at  FS  71.5 
in  the  FPG  test  section.  Figure  5.6  presents  a  comparison  between  the  computed  Mach 
number  profile  and  Miller’s  (16)  hot  film  data  obtained  at  the  same  location.  Agreement 
between  the  two  curves  is  excellent  within  the  boundary  layer  but  once  again  a  Mach 
number  mismatch  is  evident  in  the  freestream.  In  this  case,  the  computed  solution  over¬ 
predicts  the  freestream  Mach  number  at  this  flow  station  by  less  than  1.2%. 

Comparisons  between  the  computed  and  measured  turbulent  shear  stress  at  this  flow 
station  (Figure  5.7)  reveal  a  large  mismatch  between  the  numerically  and  experiment2dly 
generated  curves.  The  character  of  the  experimental  data  is  captured  in  the  computational 
curve,  but  the  average  difference  in  magnitude  is  greater  than  50%,  with  the  numerical 
results  severely  under-predicting  the  experimental  data. 

Comparison  with  the  LDV  data  of  Luker  (15)  at  this  flow  station  reveals  similar, 
though  not  identical,  results.  The  computed  non-dimensional  velocity  profile  shown  in 
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Figure  5.4  Nondimensional  Turbulent  Shear  Stress  vs  y/^„,  ZPG,  FS  60.0/71.5 


Figure  5.5  k  vs  y/^u,  ZPG,  FS  60.0/71.5 
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Figure  5.8  tracks  with  the  experimental  data  from  the  wall  to  the  freestream,  accurately 
predicting  every  aspect  of  the  profile.  Agreement  with  the  experimentally  generated  nondi- 
mensional  turbulent  shear  stress  curve,  however,  is  not  as  good.  As  shown  in  Figure  5.9, 
the  computed  curve  over-predicts  the  experimental  data  by  roughly  50%  while  maintaining 
the  character  of  the  experimentally  generated  profile.  A  comparison  of  Figures  5.9  and  5.7 
reveals  that  large  discrepancies  exist  between  the  data  sets  generated  by  Miller  (16)  and 
Luker  (15).  Luker  (15)  has  suggested  that  this  discrepancy  may  be  due  to  the  failure  of  the 
assumption  that  p'  =  0,  an  assumption  required  in  the  hot-film  data  reduction  technique. 
While  this  assumption  is  considered  valid  in  zero  pressure  gradient  regions,  Luker  (15) 
asserts  that  it  may  not  be  valid  when  a  pressure  gradient  is  present.  The  computational 
profile  splits  the  difference  between  the  hot-film  and  LDV  measurements.  Figure  5.10 
compares  the  computed  and  experimental  k  profiles  at  this  flow  station.  As  was  seen  in 
the  ZPG  test  section,  agreement  between  the  computational  and  experimental  profiles  is 
outstanding  through  the  majority  of  the  boundary  layer.  Differences  near  the  wall  may 
be  attributed  to  the  breakdown  of  an  assumption  made  by  Luker  (15)  and  Hale  (9)  in  the 
calculation  of  the  experimental  k  data  points  in  regions  of  FPG.  An  assumption  regarding 
the  magnitude  of  tn'  is  required  to  calculate  the  turbulent  kinetic  energy  because  although 
the  LDV  setup  did  not  allow  for  the  measurement  of  velocity  components  in  the  transverse 
direction,  fluctuations  in  this  direction  still  contribute  to  k.  The  assumption  made  in  the 
experimental  analysis  is  that  w'  =  v'. 

5.1.3  APG  Test  Section.  Profiles  presented  in  this  section  were  generated  at  FS 
71.0  (Dotter  (5))  and  FS  68.0  (Hale  (9)).  As  previously  mentioned,  FS  71.0  in  the  APG 
test  section  actually  lies  in  a  FPG  region.  This  should  not  affect  comparisons  between  the 
computational  results  and  those  presented  by  Dotter  (5),  although  the  profiles  will  not 
necessarily  represent  those  characteristic  of  an  APG  flowfield. 

Figure  5.11  is  a  plot  of  the  local  Mach  number  versus  y/Su  at  FS  71.0.  In  this  figure 
a  sizeable  mismatch  is  seen  between  the  computed  and  calculated  Mach  number  profiles. 
The  profile  itself  is  similar,  but  the  magnitude  of  the  computed  curve  is  appreciably  lower 
than  that  of  the  experimental  data.  The  error  appears  almost  as  an  “offset”  type  of  error 
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Figure  5.10  k  vs  y/6u,  FPG,  FS  71.5 

and  may  be  due  to  the  upstream  flowfield  specification  as  mentioned  previously.  In  the 
freestream,  the  magnitude  of  this  error  is  3.2%. 

As  may  be  seen  in  Figure  5.12,  computed  nondimensional  turbulent  shear  stress  pro¬ 
files  at  this  flow  station  are  not  in  close  agreement  with  the  hot  film  data  of  Dotter  (5), 
although  some  of  the  characteristics  of  the  experimental  data  are  present.  The  computa¬ 
tional  curve  exhibits  a  linear  region  from  0.3  <  y/6  <  0.7,  the  slope  of  which  is  only  slightly 
greater  than  that  of  the  linear  portion  of  the  hot  film  curve.  Additionally,  the  experimental 
profile  shows  a  sharp  increase  in  turbulent  shear  as  the  wall  is  approached.  This  trend  is 
also  seen  in  the  numerical  results,  the  difference  once  again  being  that  the  magnitude  of 
the  computed  curve  is  significantly  less  than  that  of  its  experimental  counterpart. 

Comparisons  between  the  numerical  curves  and  the  LDV  data  of  Hale  (9)  gener¬ 
ated  at  FS  68.0  appear  more  promising.  The  nondimensional  velocity  profiles,  shown 
in  Figure  5.13,  agree  extremely  well  through  the  boundary  layer.  The  computed  turbu¬ 
lent  shear  stress  curve  (Figure  5.14)  displays  the  same  slope  ais  the  experimental  data,  but 
under-predicts  the  magnitude  of  the  stress.  The  computational  k  profile  is  similarly  under- 


predicted,  as  may  be  seen  in  Figure  5.15,  where  the  character  of  the  experimental  profile 
is  replicated  very  well  numerically;  only  the  magnitude  is  mis-represented.  As  mentioned 
previously,  the  mismatch  seen  in  these  last  two  figures  may  be  the  result  of  difficulties  in 
determining  the  boundary  layer  thickness,  which  is  compounded  at  this  flow  station  by 
the  coalescence  of  pressure  waves  near  y  =  6  (see  Figure  5.47).  This  coalescence  has  the 
effect  of  smearing  the  edge  of  the  boundary  layer,  making  determination  of  the  boundary 
layer  thickness  a  difficult  exercise.  Altering  the  boundary  layer  thickness  utilized  in  the 
generation  of  these  curves  will  shift  the  profiles  vertically,  diminishing  both  the  apparent 
accuracy  displayed  in  Figure  5.13  and  the  discrepancies  seen  in  Figures  5.14  and  5.15. 


Figure  5.13  u/ue  vs  y/S^,  APG,  FS  68.0 


5.2  Comparison  to  Baldwin- Lomax  Results 

In  order  to  further  evaluate  the  effectiveness  of  the  k  -  w  turbulence  model  and  its 
present  implementation,  comparisons  to  numerical  solutions  obtained  through  the  use  of 
the  Baldwin- Lomax  turbulence  model  were  conducted.  This  section  presents  a  comparison 
between  solutions  obtained  with  these  two  models  for  the  three  test  geometries,  utilizing 
the  LDV  data  of  Luker  (15)  and  Hale  (9)  as  the  basis  for  comparison. 
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Figure  5.14  Nondimensional  Turbulent  Shear  Stress  vs  y/^ut  APG,  FS  68.0 


Figure  5.15  k  vs  y/^u»  APG,  FS  68.0 
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5.2.1  Zero  Pressure  Gradient  (ZPG)  Test  Section.  Figure  5.16  is  a  plot  of  the 
numerically  and  experimentally  determined  nondimensional  velocity  profiles  at  FS  60.0 
(numerical  data)  and  FS  71.5  (ZPG  LDV  data).  As  may  be  seen  in  this  figure,  agreement 
between  the  two  numerical  curves  in  the  near- wall  and  wake  regions  is  good.  The  profiles 
differ  over  a  larger  portion  of  the  boundary  layer,  extending  from  0.05  <  yjS  <  0.7.  Similar 
disagreement  was  observed  in  Section  4.5,  where  the  Baldwin- Lomax  model  under-predicts 
the  velocity  in  the  logarithmic  region  of  the  wall  law  plot.  Agreement  between  the  Baldwin- 
Lomax  and  LDV  velocity  profiles  at  this  flow  station  is  therefore  inferior  to  that  observed 
with  the  k  —  u  solution. 


iVU, 

Figure  5.16  u/u*  vs  y/6u,  ZPG,  FS  60.0/71.5 


Turbulent  shear  stress  profiles  taken  at  this  flow  station  (see  Figure  5.17)  reveal 
similar  trends.  The  Baldwin-Lomax  profiles  agree  well  with  the  k  —  u  solution  both  near 
the  wall  and  in  the  far  field,  but  outside  of  these  regions,  the  turbulent  shear  stress  values 
predicted  by  the  the  Baldwin-Lomax  model  are  marginally  higher  than  those  of  the  k  —  u 
model.  Once  again,  the  character  of  the  experimental  profile  is  better  approximated  by 
the  profile  generated  through  the  use  of  the  k  —  u  turbulence  model. 
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Figure  5.17  Nondimensional  Turbulent  Shear  Stress  vs  y/^u>  ZPG,  FS  60.0/71.5 

The  disparity  between  the  Baldwin- Lomax  and  k—u  profiles  seen  in  the  figures  above 
must  be  linked  to  the  computed  eddy  viscosity  distribution,  since  this  is  the  method  by 
which  the  effects  of  turbulent  fluctuations  are  incorporated  into  the  mean  flow  equations. 
Figure  5.18  contains  profiles  of  nrlti  vs  yjS^  for  both  the  Baldwin-Lomax  and  k~u  com¬ 
putations,  taken  at  FS  60.0.  As  may  be  seen  in  this  figure,  the  profiles  differ  substantially 
in  the  region  between  0.1  <  y|6^^  <  0.7,  where  the  Baldwin-Lomax  model  under-predicts 
the  eddy  viscosity  by  a  large  margin.  This  region  is  precisely  that  in  which  the  nondi¬ 
mensional  velocity  profiles  (see  Figure  5.16)  disagree.  Disagreement  between  the  turbulent 
shear  stress  profiles  examined  above  is  limited  to  the  region  given  by  0.1  <  y/rf„  <  0.4,  but 
the  largest  offset  (see  Figure  5.17)  occurs  at  approximately  the  same  y/6„  location  as  the 
largest  offset  in  the  eddy  viscosity  profiles. 

Notice  that  the  k  -  lj  eddy  viscosity  profile  presented  in  Figure  5.18  does  not  go  to 
0  in  the  freestream.  This  is  the  result  of  the  “experimental  profile”  upstream  boundary 
conditions,  which  specified  a  finite  eddy  viscosity  in  the  freestream.  Physically,  this  may 
be  thought  of  as  the  result  of  finite  freestream  turbulence  intensities.  Because  the  k  -  u 
turbulence  model  determines  the  eddy  viscosity  based  on  the  state  of  the  entire  flowfield. 
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this  upstream  specification  is  reflected  in  the  downstream  profiles.  As  with  all  algebraic 
turbulence  models,  the  Baldwin-Lomax  model  relies  only  on  the  local  state  of  the  flow  to 
determine  the  eddy  viscosity;  the  upstream  boundary  condition  specification  therefore  has 
little  or  no  effect  on  the  downstream  eddy  viscosity  distribution. 


Figure  5.18  firln  vs  ZPG,  FS  60.0/71.5 


5.2.2  Favorable  Pressure  Gradient  (FPG)  Test  Section.  This  section  compares 
numerical  profiles  generated  in  the  FPG  test  section  at  FS  71.5  with  the  LDV  data  of 
Luker  (15)  obtained  at  the  same  flow  station.  Figures  5.19  and  5.20  are  plots  of  the  nondi- 
mensional  velocity  and  turbulent  shear  stress,  respectively,  at  this  location.  These  figures 
exhibit  the  same  trends  with  regards  to  the  Baldwin-Lomax  solutions  that  were  observed 
in  the  ZPG  region,  with  under-predictions  of  the  velocity  and  over-predictions  in  the  tur¬ 
bulent  shear  stress  appearing  between  0.05  <  y/6„  <  0.7.  As  before,  the  k  -  u  turbulence 
model  appears  to  more  accurately  represent  the  experimental  data.  Examination  of  the 
eddy  viscosity  profiles  at  this  flow  station  (see  Figure  5.21)  reveals  significant  disparities 
throughout  the  boundary  layer.  While  the  peak  eddy  viscosity  is  similarly  predicted  by 
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both  models,  the  Baldwin- Lomajc  model  predicts  that  this  peak  occurs  near  the  boundary 
layer  edge,  whereas  the  k  —  u  model  predicts  its  location  at  «  0.6. 


Figure  5.19  u/ti*  vs  jz/^ui  FPG,  FS  71.5 


5.2.3  Adverse  Pressure  Gradient  (APG)  Test  Section.  In  this  section,  computed 
nondimensional  velocity  and  turbulent  shear  stress  profiles  at  FS  68.0  in  the  APG  test 
section  are  compared;  the  LDV  data  of  Hale  (9)  serves  as  the  basis  for  comparison.  In  Fig¬ 
ure  5.22,  the  Baldwin- Lomax  solution  is  seen  to  deviate  only  slightly  from  the  k  —u  profile. 
Unfortunately,  this  profile  deviation  is  such  that  the  computed  curve  no  longer  matches 
the  character  of  the  experimental  data.  Similarly,  in  Figure  5.23  the  overall  magnitude 
of  the  turbulent  shear  stress  is  similar  to  that  of  the  k  —  u  solution,  but  the  character  of 
the  experimental  curve  is  lost.  Inspection  of  Figure  5.24  reveals  good  agreement  between 
the  computed  eddy  viscosity  profiles  in  the  inner  and  outer  regions  of  the  boundary  layer, 
but  a  significant  under- prediction  of  this  quantity  between  0.05  <  y/6„  <  0.6.  This  is 
seen  as  the  root  of  the  large  disparities  in  the  velocity  and  turbulent  shear  stress  profiles 
mentioned  previously. 
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Figure  5.20  Nondimensional  Turbulent  Shear  Stress  vs  y/^u>  FPG,  FS  71.5 


Figure  5.21  vs  y,  FPG,  FS  71.5 


Figure  5.22  «/ti<  vs  y/^u,  APG,  FS  68,0 


Figure  5.23  Nondimensional  Turbulent  Shear  Stress  vs  y/^u,  APG,  FS  68.0 
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Figure  5.24  tirln  vs  y,  APG,  FS  68.0 


5.2.4  Surface  Effects.  Figure  5.25  displays  the  distribution  of  the  wall  shear 
stress  along  the  test  sections  for  the  Baldwin-Lomax  and  k  —  u  computations.  Neglecting 
the  entrance  region,  the  Baldwin-Lomax  turbulence  model  predicts  wall  shear  values  in 
the  ZPG  region  that  are  8%  to  13%  less  than  those  predicted  using  the  k  —  u  turbulence 
model.  Using  the  correlation  from  van  Driest  II  theory  (24)  (1),  Luker  (15)  calculated 
the  wall  shear  stress  at  FS  71.5  in  the  ZPG  test  section  as  67.0  Pa.  Extrapolating  the 
slope  of  the  r^aii  vs  x  curve  from  upstream  of  the  test  sections  to  71.5  cm  shows  that  the 
computed  r^aii  at  this  station  would  be  approximately  66.8  Pa,  a  difference  of  0.3%  from 
the  experimental  value.  A  similar  extrapolation  performed  on  the  Baldwin-Lomax  curve 
yields  a  wall  shear  stress  of  64  Pa.  Interestingly,  the  Baldwin-Lomax  results  indicate  that 
>  0,  a  trend  that  predicted  by  neither  the  k  —  u  results  nor  theory.  Predictions  of 
waU  shear  in  the  FPG  test  section  are  similar  in  both  slope  and  magnitude,  whereas  the 
Baldwin-Lomax  computations  in  the  APG  test  section  under-predict  the  wall  shear  by  a 
large  margin  when  compared  to  the  k  —  u  solution. 

Figure  5.26  contains  plots  of  the  computed  surface  pressure  as  a  function  of  location. 
The  pressure  distribution  predicted  through  the  use  of  the  Baldwin-Lomax  model  is  very 
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similar  to  that  obtained  using  the  Ar  -  w  model  over  the  majority  of  the  tunnel  surface 
for  both  experimental  test  sections.  This  agreement  breaks  down  somewhat  in  the  region 
between  FS  58.0  and  FS  65.0.  In  this  region,  a  decrease  in  surface  pressure  is  predicted 
by  both  models,  the  artifact  of  a  series  of  pressure  waves  that  are  present  in  the  computa¬ 
tional  domain  as  the  result  of  the  boundary  conditions  applied  at  the  inflow  station.  The 
magnitude  and  duration  of  this  pressure  fluctuation  is  larger  with  the  k  -  u  model  than 
with  the  Baldwin-Lomax  model.  In  the  APG  and  FPG  test  sections,  however,  Pu,oH  and 
are  seen  to  be  elfectively  independent  of  the  turbulence  model  utilized. 


5.2.5  Computational  Expense  and  Convergence  Characteristics.  Although  the 
solutions  obtained  with  the  Baldwin-Lomax  turbulence  model  are  in  general  found  to  be 
less  accurate  than  that  obtained  through  use  of  the  k  —  u  turbulence  model,  use  of  the 
former  for  flowfield  analysis  may  be  justified  if  a  significant  savings  in  computational  time 
may  be  realized  by  doing  so.  A  comparison  of  required  Central  Processing  Unit  (CPU) 
time  for  a  fixed  number  of  iterations  using  the  two  turbulence  models  has  been  performed. 
For  this  comparison,  the  code  was  run  for  10,000  iterations  on  grid  APG2,  a  151x151  grid 


Figure  5.26  vs  x 

having  22801  cell  centers  on  a  Cray  Research  C916.  In  order  to  minimize  the  effect  of 
read  and  write  time  requirements,  no  data  output  subroutines  were  called  during  these 
computations.  The  k  —u  computations  required  3012.529  seconds  of  CPU  time  for  a  total 
of  1.3212x10"®  seconds  per  iteration  per  cell  center.  The  Baldwin-Lomax  computations 
required  2811.601  seconds,  or  1.2331x10"®  seconds  per  iteration  per  cell  center,  a  6.67% 
reduction  in  CPU  usage  over  the  k  —  u  model. 

The  convergence  histories  of  computations  utilizing  the  Baldwin-Lomax  model  are 
significantly  different  than  those  exhibited  when  k  —  u  modeling  is  used.  Figure  5.27 
shows  characteristic  convergence  histories  for  Baldwin-Lomax  and  k  -  tv  computations 
performed  on  Grid  FPG2.  The  CFL  number  ^used  for  these  calculations  was  0.5.  The 
convergence  history  typical  of  the  k  —  u  turbulence  model  displays  an  initial  increase  in 
the  residual,  followed  by  a  long  steady  decrease  of  roughly  12  orders  of  magnitude  before 
bottoming  out  near  2.5x10"“.  Convergence  to  this  level  is  seldom  required,  however, 
since  the  solutions  effectively  cease  to  evolve  after  the  residual  drops  four  or  five  orders 
of  magnitude.  The  convergence  history  characteristic  of  the  Baldwin-Lomax  turbulence 
model  shows  a  much  smaller  initial  jump  in  the  residual  followed  by  a  shallow  decrease 
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of  roughly  4  orders  of  magnitude.  At  this  point  oscillations  in  the  solution  preclude  any 
further  reduction  in  the  residual  unless  the  turbulence  model  is  frozen,  as  shown  in  this 
figure.  These  oscillations  are  not  eliminated  by  a  reduction  in  the  CFL  number  and  are 
thought  to  be  linked  to  the  specification  of  and  within  the  Baldwin-Lomax 
subroutine  (see  Section  3.3.1).  As  ymax  bounces  between  two  adjacent  nodes,  the  flowfield 
is  continually  attempting  to  adjust  to  the  new  formulation  for  fir  at  that  plane  within 
the  flow.  This  behavior  may  be  at  least  partially  eliminated  by  the  introduction  of  finer 
computational  meshes,  but  the  intent  of  this  section  is  to  evaluate  the  performance  of 
the  models  on  the  same  mesh  for  the  same  flowfield  conditions.  For  the  purposes  of  this 
investigation,  Baldwin-Lomax  computations  are  deemed  converged  when  this  oscillation 
continues  for  roughly  10000  iterations  unchanged.  Fortunately,  this  oscillatory  behavior 
appears  after  the  flowfield  properties  of  interest  {T^aihPwaiuV’iT^y)  have  ceased  to  change. 
It  should  be  noted  that  this  oscillation,  though  typical  of  Baldwin-Lomax  computations, 
did  not  appear  in  every  such  computation.  Variations  in  freestream/upstream  conditions 
had  a  profound  effect  on  the  convergence  history  of  Baldwin-Lomax  computations;  in  fact, 
several  Baldwin-Lomax  computations  (using  different  flow  conditions)  converged  to  the 
same  degree  as  the  companion  k  —  u  computations. 

In  summary,  the  Baldwin-Lomax  turbulence  model  produces  results  for  these  flow- 
fields  that  are  found  to  be  less  accurate  and  more  sensitive  to  the  inflow  conditions  than 
those  obtained  through  use  of  the  k  —  u)  turbulence  model.  This  may  be  the  result  of  the 
Baldwin-Lomax  model’s  inability  to  account  for  flow  history  effects,  an  ability  possessed 
by  the  k  —  u  model,  as  mentioned  in  Chapter  1.  As  such,  use  of  this  model  based  on 
potential  CPU  time  savings  is  not  a  viable  option  for  flowfields  of  the  type  examined  in 
this  study.  As  an  alternative,  implicit  solution  of  the  k  —  u  equations  has  the  potential  to 
significantly  reduce  computational  expense  while  maintaining  the  high  level  of  accuracy  of 
the  k  —  ijj  formulation. 

5.3  Evaluation  of  Inflow  Boundary  Condition  Effects 

The  uniqueness  of  the  solution  of  a  system  of  differential  equations  is  guaranteed  by 
the  form  of  the  boundary  conditions  imposed  on  that  system  of  equations.  This  mathe- 
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Herations 

Figure  5.27  Convergence  History,  k  —  u  vs  Baldwin' Lomax 

matical  reality  drives  the  need  to  evaluate  the  effect  of  boundary  condition  variation  on 
the  solutions  generated  in  this  effort.  Since  the  flowfield  being  studied  is  for  the  most  part 
supersonic,  the  preferred  direction  of  information  propagation  is  largely  in  the  streamwise 
direction.  Thus,  variations  in  the  upstream  boundary  conditions  may  be  expected  to  have 
a  significant  effect  on  the  flowfield. 

The  two  sets  of  boundary  conditions  used  in  this  investigation  are  described  in  Sec¬ 
tion  4.6.  The  following  sections  compare  and  contrast  the  solutions  obtained  with  each  set 
of  boundary  conditions,  once  again  using  the  LDV  data  of  Luker  (15)  and  Hale  (9)  as  the 
basis  for  comparison. 

5.S.1  Zero  Pressure  Gradient  (ZPG)  Test  Section.  This  section  compares  the 
numerical  solutions  obtained  at  FS  60.0  with  the  LDV  data  extracted  at  FS  71.5  in  the 
ZPG  test  section.  Figure  5.28  is  a  plot  of  the  nondimensional  velocity  as  a  function  of 
y/6  through  the  boundary  layer  at  this  flow  station.  Although  the  numerically  generated 
solutions  agree  closely  in  the  wall  region,  this  agreement  disappears  above  y/6  «  0.2 
Examination  of  the  freestream  inflow  curve  reveals  a  velocity  deficit  in  the  upper  two- 
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thirds  of  the  boundary  layer  that  is  not  seen  in  either  the  experimental  profile  inflow  curve 
or  the  experimental  data.  This  velocity  deficit  may  be  linked  to  an  over-prediction  of  the 
turbulent  shear  stress  in  the  same  region,  as  is  seen  in  Figure  5.29. 


Figure  5.28  u/ue  vs  ZPG,  FS  60.0/71.5 


5.3.2  Favorable  Pressure  Gradient  (FPG)  Test  Section.  Profiles  in  this  section 
were  generated  at  FS  71.5  in  the  FPG  test  section.  Figure  5.30  shows  that  the  nondi- 
mensionaJ  computed  velocity  profiles  differ  only  slightly  at  this  flow  station;  both  do  an 
acceptable  job  approximating  the  experimental  data.  Similarly,  Figure  5.31  shows  close 
agreement  between  the  computed  solutions.  Comparing  the  features  of  Figures  5.30  and 
5.31  reveals  that  for  a  given  yjb  location,  high  nondimensional  velocity  is  typically  coupled 
with  reduced  turbulent  shear  stress.  Thus,  an  over-prediction  of  the  local  velocity  appears 
to  be  linked  to  an  under-prediction  of  the  turbulent  shear  stress.  This  trend  is  also  present 
in  the  ZPG  test  section,  as  may  be  seen  in  Figures  5.28  and  5.29. 

5.3.3  Adverse  Pressure  Gradient  (APG)  Test  Section.  Profiles  presented  in 
this  section  were  generated  at  FS  68.0  in  the  APG  test  section.  Agreement  between  the 
numerical  profiles  is  similar  to  that  which  was  seen  in  the  ZPG  test  section.  The  differences 
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Figure  5.31  Nondimensional  Turbulent  Shear  Stress  vs  FPG,  FS  71.5 

may  once  again  be  characterized  by  a  velocity  deficit  and  accompanying  turbulent  shear 
stress  surplus  between  0.2  <  yl6  <  0.75.  The  velocity  profile  is  once  again  more  accurately 
predicted  numerically  when  the  experimental  profile  boundary  conditions  are  used.  The 
turbulent  shear  stress  profile  generated  through  the  use  of  the  freestream  inflow  conditions, 
however,  more  accurately  predicts  the  magnitude  of  this  quantity  over  a  fairly  large  region 
of  the  boundary  layer.  Unfortunately,  this  agreement  is  restricted  to  magnitude  effects  only; 
the  general  shape  of  the  turbulent  shear  stress  profile  is  still  more  accurately  reproduced  by 
the  computations  utilizing  the  experimental  profile  conditions.  Finally,  the  converse  of  the 
trend  observed  in  the  above  sections  appears  to  hold  as  well;  numerical  under-prediction 
of  the  velocity  is  coupled  with  an  over-prediction  of  the  turbulent  shear  stress. 

5.3.4  Surface  Effects.  The  effect  of  the  upstream  boundary  condition  variation  is 
not  limited  to  the  velocity  and  shear  stress  profiles.  Differences  exist  in  the  surface  pressure 
and  wall  shear  stress  distributions  as  well.  Figure  5.34  shows  that  while  the  computed  w^l 
shear  stress  values  are  the  same  for  both  models  in  the  upstream  ZPG  region,  wall  shear 
predictions  in  the  pressure  gradient  regions  differ.  This  figure  indicates  that  the  magnitude 
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of  the  change  in  the  wall  shear  stress  in  the  pressure  gradient  regions  is  smaller  with  the 
freestream  inflow  conditions  than  that  obtained  through  the  use  of  the  experimental  profile 
inflow  conditions. 


I 


Figure  5.34  r^aii  vs  x 


Figure  5.35  shows  that  while  the  pressure  gradient  within  the  test  sections  is  similarly 
predicted  by  both  inflow  boundary  conditions,  a  fairly  large  discrepancy  exists  in  the 
computed  surface  static  pressure  distributions.  In  this  figure,  computations  using  the 
freestream  inflow  conditions  predict  significantly  higher  static  pressures  thsm  those  utilizing 
the  experimental  inflow  profile.  Interestingly,  it  is  the  curve  generated  using  the  freestream 
inflow  conditions  that  most  closely  approximates  the  experimental  data  (Luker  (15),  FPG 
data  and  Hale  (9),  APG  data).  Recall  that  both  inflow  conditions  use  essentially  the  same 
data,  the  difference  being  that  the  latter  utilizes  more  flowfield  information  by  specifying 
nonuniform  velocity,  temperature,  density,  k,  and  u;  profiles.  Thus,  there  is  no  reason  to 
expect  that  the  pressure  distribution  should  vary  when  one  set  of  boundary  conditions  is 
used  in  lieu  of  the  other. 

The  differences  observed  may  be  partially  explained  through  an  examination  of  Fig¬ 
ures  5.36  and  5.37.  The  lower  frame  of  these  figures  contains  a  numerical  contour  plot  of 
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Figure  5.35  vs  * 

the  static  pressure  in  the  test  section;  the  upper  frame  is  a  line  plot  of  the  static  pressure 
at  the  surface.  Comparing  these  figures  with  Figure  5.35  reveals  that  while  the  pressure 
at  the  inflow  plane  is  effectively  the  same  for  both  sets  of  boundary  conditions,  a  series  of 
reflected  pressure  waves  present  in  the  numerical  solution  increases  the  surface  (and  mean 
flow)  pressure  substantially  in  the  streamwise  direction.  These  pressure  waves  appear  to 
originate  at  the  leading  edge  of  the  computational  mesh  when  the  previously  undisturbed 
freestream  is  subjected  to  no-slip  boundary  conditions  at  the  inlet.  This  is  an  interesting 
feature,  but  it  does  not  explain  the  agreement  between  the  freestream  inflow  case  and  the 
experimental  data. 

Looking  to  the  experimental  setup  for  potential  answers,  Schlieren  images  taken  in 
the  test  sections  reveal  the  presence  of  “seam  shocks,”  oblique  shock  waves  which  originate 
upstream  of  the  measurement  stations  due  to  misalignment  of  the  wind  tunnel  sections. 
The  presence  of  these  weak  shock  waves  raises  two  issues:  first,  it  challenges  the  assumption 
that  the  total  pressure  is  constant  from  the  settling  chamber  to  FS  44.0,  the  location  at 
which  the  upstream  boundary  conditions  were  measured.  The  effect  of  these  seam  shocks 
would  be  to  reduce  the  total  pressure  at  that  station,  leading  to  further  reductions  in  the 
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static  pressure.  Thus,  the  static  pressure  assumed  at  FS  44.0  is  the  maximum  possible  static 
pressure  at  that  flow  station  given  the  velocity  profile  and  total  pressure  and  temperature 
measured  in  the  settling  chamber.  The  second  is  that  these  shocks  may  be  responsible 
for  raising  the  static  pressure  from  the  “assumed”  upstream  value  of  roughly  7200  Pa  to 
the  8400  Pa  measurement  obtained  by  Hale  (9)  at  FS  63.0.  In  that  case,  the  agreement 
between  the  freestream  inflow  calculation  and  the  experimental  data  is  fortuitous  at  best; 
the  non-physical  numerical  pressure  fluctuations  simply  combine  in  such  a  manner  as  to 
mimic  the  pressure  rise  in  the  actual  flow  due  to  these  unintended,  but  physical,  oblique 
shocks. 

These  reflected  pressure  waves  are  also  present  to  some  extent  in  Figures  5.38  and 
5.39,  which  were  generated  through  the  use  of  the  experimental  profile  inflow  boundary 
conditions.  The  appearance  of  pressure  fluctuations  in  these  solutions  tends  to  dispel  the 
notion  that  the  instantaneous  application  of  no-slip  at  the  leading  edge  is  responsible  for 
the  computational  pressure  fluctuations.  A  potential  link  between  these  two  calculations 
is  that  in  both  cases  the  v  component  of  velocity  is  set  to  0.  For  the  freestream  inflow 
case,  this  is  a  good  approximation,  consistent  with  the  assumption  of  an  undisturbed 
freestream.  In  the  experimental  profile  case,  however,  this  assumption  is  fairly  inaccurate; 
boundary  layer  flows  typically  have  non-zero  v  components  of  velocity.  It  is  possible  that 
this  inconsistency  is  at  least  partially  responsible  for  the  appearance  of  the  pressure  waves 
in  that  calculation. 

Finally,  a  portion  of  the  observed  pressure  mismatch  may  be  linked  to  the  stagnation 
temperature  and  pressure  utilized  in  the  specification  of  both  sets  of  boundary  conditions 
(see  Section  4.6.  For  example,  a  1  degree  difference  in  the  total  temperature  specification 
leads  to  a  140  Pa  oflfset  in  the  static  pressure.  Similarly,  uncertainty  in  the  stagnation 
pressure  measurements  reported  by  Luker  (15)  lead  to  an  88  Pa  window  over  which  the 
static  pressure  may  vary. 

5.3.5  Summary.  Specification  of  the  upstream  boundary  condition  appears  to 
have  a  significant  effect  on  the  velocity  and  shear  stress  profiles  in  the  ZPG  and  APG 
test  sections.  This  effect  is  diminished  appreciably  in  the  FPG  test  section.  Thus,  the 
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stabilizing  effect  of  a  FPG  tends  to  reduce  the  effects  of  variations  in  upstream  boundary 
condition  specification,  whereas  the  effect  of  an  APG  is  to  maintain,  if  not  amplify,  these 
differences.  From  a  practical  standpoint,  this  observation  implies  that  extra  care  must 
be  taken  when  specifying  the  upstream  boundary  conditions  for  a  flowfield  characterized 
by  ZPG  or  an  APG.  The  converse  appears  to  hold  as  well  for  FPG  flowfields.  Finally, 
accurate  specification  of  the  upstream  profile,  to  include  thermophysical  quantities  and 
all  components  of  the  velocity  vector,  is  required  to  obtain  accurate  absolute  pressure 
distributions,  although  the  magnitudes  of  pressure  gradients  are  reproduced  accurately 
regardless  of  the  inflow  condition  used. 

5.4  Flowfield  Analysis 

This  section  presents  the  results  of  a  flowfield  analysis  performed  using  the  exper¬ 
imental  profile  inflow  conditions  and  the  k  —  u  turbulence  model,  the  combination  of 
turbulence  model  and  upstream  boundary  conditions  found  to  provide  the  most  accurate 
overall  solutions.  As  in  the  previous  sections,  this  section  focuses  on  the  effects  of  the 
favorable  and  adverse  pressure  gradients  through  an  examination  of  the  velocity,  turbulent 
shear  stress,  and  eddy  viscosity  profiles,  the  static  pressure  distribution  along  the  surface, 
and  the  distribution  of  Ty,an  through  the  test  sections. 

5.4-1  Favorable  Pressure  Gradient  (FPG)  Test  Section.  The  flowfield  in  the  FPG 
test  section  is  characterized  by  a  progressive  and  mild  expansion  beginning  at  FS  65.0  and 
continuing  to  the  end  of  the  computational  domain.  A  numerical  Schlieren  image  of  this 
test  section,  showing  density  gradients  in  the  y  direction,  |^,  is  presented  in  Figure  5.40(a). 
In  this  figure,  flow  is  from  left  to  right  and  lightly  shaded  areas  represent  regions  of  positive 
density  gradients;  dark  areas  correspond  to  negative  gradients.  The  expansion  is  visible 
as  an  incrementally  brighter  triangular  region  in  the  lower  right  hand  comer  of  the  figure. 
In  this  figure,  the  boundary  layer  at  the  lower  surface  is  seen  to  respond  immediately 
to  the  presence  of  this  pressure  gradient  by  increasing  in  thickness,  while  the  boundary 
layer  at  the  upper  surface  appears  to  be  unaffected  by  the  expansion.  This  observation  is 
consistent  with  those  of  Spina,  Smits  &  Robinson  (21),  who  found  that  “...bulk  dilitation 
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decreases  the  wall  shear  stress  and  increases  the  boundary-layer  thickness.”  Referring  back 
to  Figure  5.38,  the  average  pressure  gradient  at  the  lower  surface  in  the  expansion  region 
is  fairly  constant  and  may  be  calculated  as  —4.5728x10'*  Pa/m. 


(a)  (b)  Luker  (15) 

Figure  5.40  Schlieren  Imagery  of  FPG  Test  Section 


Figure  5.40(b)  is  an  actual  Schlieren  photograph  of  the  FPG  test  section.  Pronainent 
in  this  figure  are  the  “seam  shocks”  mentioned  in  Section  5.3.4.  One  appears  as  a  light 
line  passing  through  the  frame  from  upper  left  to  lower  right,  the  other  begins  at  the  lower 
left,  moving  towards  the  upper  right  corner  but  dissipating  away  before  it  arrives  there. 
This  pattern  is  initially  symmetric  about  the  centerline  of  the  tunnel.  In  this  figure,  a 
thickening  of  the  boundary  layer  is  evident  as  the  flow  passes  into  the  expansion  region, 
but  the  expansion  is  more  difficult  to  visualize  due  to  noise  in  the  flow  and  the  the  arrival 
of  the  first  seam  shock. 

The  remainder  of  this  section  will  concentrate  on  profiles  obtained  at  different  flow 
stations  within  the  FPG  test  section.  The  location  of  these  stations  and  the  orientation 
of  the  cutting  planes  (lines)  with  respect  to  the  test  section  are  shown  in  Figure  5.41,  a 
composite  figure  based  on  a  contour  plot  of  the  pressure  gradient  in  the  x  direction, 
Superimposed  over  this  contour  plot  are  the  planes  along  which  data  has  been  extracted; 
directly  above  the  plot  is  a  representation  of  the  surface  pressure  distribution  at  the  bottom 
wall  in  the  test  section.  Immediately  evident  in  this  contour  plot  are  the  reflected  pressure 
waves  mentioned  previously  and  a  significantly  more  visible  expansion  region. 
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Figure  5.42  Velocity  vs  y,  FPG  Test  Section 


Figure  5.43  Nondimensional  Turbulent  Shear  Stress  vs  FPG  Test  Section 
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of  the  profiles  as  the  cutting  planes  are  moved  downstream.  This  figure  shows  that  while 
neither  the  peak  turbulent  kinetic  energy  (TKE)  near  the  wall  nor  the  area  under  the 
k  -  y  curve  are  significantly  reduced,  the  distribution  of  TKE  through  the  boundary  layer 
is  altered  by  the  presence  of  the  FPG.  This  alteration  consists  of  a  shift  which  increases 
the  TKE  near  the  edge  of  the  boundary  layer  while  reducing  its  magnitude  close  to  the 
wall.  A  physical  argument  is  proposed:  High  inertia  freestream  flow  resists  the  influence 
of  the  curvature  of  the  wall,  increasing  the  boundary  layer  thickness.  As  a  result  of  this 
increased  boundary  layer  thickness,  the  velocity  gradient  near  the  wall  decreases,  resulting 
in  decreased  shear  stress  at  the  wall  and  decreased  energy  transfer  to  the  wall  through 
viscous  effects.  This  reduction  in  the  wall  shear  stress  may  be  clearly  seen  in  Figure  5.34, 
where  the  shear  stress  at  the  wall  is  reduced  from  72  Pa  to  46  Pa  over  the  course  of  8  cm. 
Since  the  total  energy  in  the  flow  is  conserved,  the  reduction  in  energy  transfer  to  the  wall 
must  be  balanced  by  an  increase  in  the  energy  transfer  elsewhere.  This  shifts  the  TKE 
curve  toward  the  edge  of  the  boundary  layer,  where  increased  transfer  of  TKE  into  the 
mean  flow  accounts  for  the  reduction  in  energy  transferred  to  the  wall. 


k  [m’/s*] 


Figure  5.44  Turbulent  Kinetic  Energy  vs  y,  FPG  Test  Section 
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Figure  5.45  shows  the  evolution  of  the  eddy  viscosity  profiles  through  the  FPG  test 
section.  The  peak  eddy  viscosity  is  seen  to  decrease  and  move  away  from  the  wall  as  the 
cutting  plane  is  moved  downstream.  This  magnitude  reduction  and  peak  shift  are  thought 
to  be  responsible  for  the  reduction  in  the  turbulent  shear  stress  seen  in  Figure  5.43  and 
the  re-distribution  of  turbulent  kinetic  energy  seen  in  Figure  5.44.  Examination  of  the 
eddy  viscosity  levels  in  the  freestream  reveals  that  these  values  are  also  reduced  under 
the  influence  of  the  FPG.  As  noted  in  Spina,  Smits  &  Robinson  (21),  a  FPG  tends  to 
decrease  the  magnitude  of  turbulent  fluctuations  in  compressible  boundary  layer  flow, 
even  re-laminarizing  these  flowfields  if  the  FPG  is  strong  enough.  Computationally,  this 
result  arises  from  a  reduction  in  the  eddy  viscosity,  a  feature  that  is  clearly  predicted  by 
the  k  —  u  turbulence  model. 


Figure  5.45  Hr  In  vs  y/^u>  k  —  oj  Modeling,  FPG  Test  Section 

5.4.2  Adverse  Pressure  Gradient  (APG)  Test  Section.  The  flowfield  in  the  APG 
Test  Section  displays  a  mixed  set  of  features.  As  mentioned  in  Chapter  II,  regions  of 
both  favorable  and  adverse  pressure  gradient  are  present  in  this  test  section.  In  the  region 
upstream  of  FS  68.0,  (see  Figure  5.39)  the  lower  wall  is  concave  up,  creating  a  region 
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of  adverse  pressure  gradient.  The  magnitude  of  the  pressure  gradient  in  this  region  is 
relatively  constant  and  calculated  to  be  1.3425x10®  Pa/m.  Beyond  FS  68.0,  the  radius 
of  curvature  increases  until  an  inflection  occurs,  transitioning  the  pressure  gradient  from 
adverse  to  favorable.  Hence,  aft  of  FS  69.0  the  pressure  gradient  is  favorable;  beyond  FS 
71.0  the  magnitude  of  this  gradient  is  fairly  constant  and  equal  to  —1.6201x10®  Pa/m 


(a)  (b)  Hale  (9) 

Figure  5.46  Schlieren  Imagery  of  APG  Test  Section 


Figure  5.46(a)  is  a  numerical  Schlieren  image  of  this  test  section.  Strong  density 
gradients  are  evident  in  the  shock  wave,  which  strengthens  near  the  end  of  the  test  sec¬ 
tion.  The  boundary  layer  is  seen  to  decrease  in  thickness  initially  behind  the  coalescence 
of  pressure  waves  near  the  wall,  and  then  to  increase  in  the  region  of  FPG  which  follows. 
Figure  5.46(b)  is  an  actual  Schlieren  image  of  this  test  section,  displaying  the  same  charac¬ 
teristics  as  the  numerical  prediction.  Note  that  the  “seam  shocks”  seen  in  Figure  5.40(b) 
are  present  in  this  image  as  well.  The  effect  of  these  weak  oblique  waves  is  diminished  in 
this  test  section,  however,  due  to  the  strength  of  the  main  shock  structure. 

Figure  5.47  is  a  contour  plot  of  the  pressure  gradient  in  the  x  direction,  over  which 
are  superimposed  the  locations  of  the  data  planes  considered  in  this  section.  From  this 
figure,  the  development  of  the  shock  wave  is  readily  visible,  transitioning  from  a  simple 
coalescence  of  pressure  waves  into  a  proper  shock  near  the  end  of  the  test  section.  As  was 
done  previously,  a  representation  of  the  surface  static  pressure  in  the  test  section  is  shown 
above  this  contour  plot. 
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Figure  5.47  Location  of  Measurement  Stations,  APG  Test  Section 


As  may  be  expected  and  seen  in  Figure  5.48,  the  combination  of  the  pressure  gradients 
and  the  shock  wave  has  a  dramatic  effect  on  the  velocity  profiles  in  the  test  section.  Initially, 
the  effect  of  the  APG  on  the  velocity  profiles  is  to  decrease  the  edge  velocity  and  thin  the 
boundary  layer,  as  may  be  seen  by  comparing  the  profiles  at  FS  66.0  and  67.0.  By  FS  68.0, 
the  shock  is  strong  enough  to  have  a  smearing  effect  on  the  boundary  layer  edge,  making 
that  edge  ill-defined.  Aft  of  FS  68.0,  the  effect  of  the  shock  is  restricted  primarily  to  the 
flow  outside  the  boundary  layer.  It  is  in  this  region  (aft  of  FS  68.0)  that  the  FPG  begins 
to  dominate  the  flowfield,  manifesting  its  presence  through  an  increase  in  the  boundary 
layer  thickness. 

The  nondimensional  turbulent  shear  stress  profiles  in  this  section  also  reflect  the 
varying  character  of  the  flowfield.  The  peak  turbulent  shear  stress  near  the  wall  is  seen  to 
increase  with  streamwise  position  initially,  maintaining  a  profile  that  resembles  that  seen 
in  the  ZPG  regions.  After  the  peak  value  is  reached  in  the  vicinity  of  FS  68.0,  the  profiles 
begin  to  change  in  character,  moving  towards  the  form  presented  in  Figure  5.43.  As  wjis 
seen  in  this  figure,  progression  along  the  FPG  region  is  accompanied  by  a  reduction  in  the 
turbulent  shear  stress  in  the  flowfield.  It  is  interesting  to  note  that  the  profiles  aft  of  FS 
68.0  in  this  test  section  do  not  immediately  take  the  form  of  those  generated  in  the  FPG 
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u(m/s] 

Figure  5.48  Velocity  vs  y,  APG  Test  Section 

test  section.  Their  shape  instead  reflects  both  APG  and  FPG  qualities.  This  observation 
serves  to  emphasize  the  importance  of  “historical”  flowfield  effects  in  turbulence  modeling. 
This  trend  is  also  reflected  in  Figure  5.34,  which  indicates  an  initial  increase  in  the  wall 
shear  stress,  followed  by  a  decrease  in  this  quantity.  Note  that  the  peak  wall  shear  occurs 
after  the  inflection  takes  place,  once  again  indicating  the  importance  of  upstream  effects. 

The  turbulent  kinetic  energy  (TKE)  profiles  shown  in  Figure  5.50  exhibit  similar 
features.  Within  the  APG  region,  the  peak  value  of  the  TKE  is  seen  to  increase  and  shift 
slightly  away  from  the  wall.  Transitioning  into  the  FPG  region  of  this  test  section,  the 
peak  value  decreases  more  dramatically,  and  the  curve  appears  to  begin  to  turn  in  on  itself 
in  an  attempt  to  return  to  the  type  of  profile  seen  earlier  in  the  FPG  test  section  (see 
Figure  5.44).  Again,  the  condition  of  the  upstream  flow  appears  to  play  a  major  role  in 
specifying  the  TKE  profiles  at  each  flow  station.  In  other  words,  the  profile  in  the  FPG 
portion  of  this  test  section  still  retains  many  of  the  features  of  the  APG  profiles,  while 
obviously  moving  towards  the  overall  shape  of  the  profiles  seen  in  the  FPG  test  section. 
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Figure  5.49  Nondimensional  Turbulent  Shear  Stress  vs  y,  APG  Test  Section 


Figure  5.50  Turbulent  Kinetic  Energy  vs  y,  APG  Test  Section 
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The  eddy  viscosity  profiles  (see  Figure  5.51)  in  the  APG  test  section  behave  in  a 
similar  manner.  From  FS  66.0  to  FS  69.0,  the  peak  eddy  viscosity  increases  at  a  fairly 
constant  rate  and  does  not  shift  vertically  within  the  boundary  layer.  The  peak  eddy 
viscosity  at  FS  70.0  is  greater  than  that  at  FS  69.0  but  has  begun  to  move  away  from 
the  wall,  as  was  seen  in  the  FPG  test  section.  This  seems  to  indicate  that  the  maximum 
eddy  viscosity  lies  somewhere  between  FS  69.0  and  70.0,  on  the  same  y  plane  as  is  seen 
between  FS  66.0  and  FS  69.0.  Then,  as  the  adverse  pressure  gradient  diminishes  in  strength 
and  the  gradient  becomes  favorable,  the  peak  eddy  viscosity  decreases  and  moves  towards 
the  boundary  layer  edge,  as  seen  at  FS  70.0  and  71.0.  As  before,  this  appears  to  be  the 
mechanism  responsible  for  shifting  the  TKE  profile  towards  the  boundary  layer  edge.  Note 
that  the  eddy  viscosity  in  the  region  between  the  shock  wave  and  the  lower  test  section 
wall  is  higher  than  that  above  the  shock,  indicating  perhaps  that  the  shock  wave  itself 
limits  the  extent  of  the  APG  to  the  region  close  to  the  wall.  In  the  region  below  the  shock 
wave  but  outside  the  boundary  layer,  eddy  viscosity  is  seen  to  increase  under  the  influence 
of  the  APG. 


Figure  5.51  ^//i  vs  y,  1:  -  u;  Modeling 
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VI.  Conclusions  and  Recommendations  for  Further  Investigation 
6.1  Conclusions 

1.  The  k  -u)  turbulence  model  provides  acceptable  computational  results  for  this  family 
of  flowfields.  Agreement  between  computed  and  experimental  velocity  and  turbu¬ 
lent  shear  stress  profiles  is  excellent,  as  is  the  agreement  between  computed  and 
experimental  wall  shear  stress  in  the  zero  pressure  gradient  test  section.  Wall  shear 
predictions  in  the  pressure  gradient  sections  appear  to  follow  published  observations, 
but  the  lack  of  experimental  data  in  these  regions  precludes  direct  comparison.  AH 
pertinent  flowfield  features  have  been  accurately  predicted  numerically  using  this 
model.  Finally,  k  —  u  computations  produced  significantly  better  results  than  the 
corresponding  Baldwin- Lomax  computations. 

2.  Flowfield  history  plays  a  significant  role  in  the  determination  of  flowfield  properties 
at  every  station  in  the  flowfield.  Velocity  and  shear  stress  profiles  evolve  under  the  in¬ 
fluence  of  both  local  and  non-local  flow  features.  Failing  to  account  for  the  non-local 
flow  features  will  introduce  errors  into  the  solution.  Although  FPG  flowfield  compu¬ 
tations  appear  to  be  fairly  insensitive  to  variations  in  upstream  boundary  conditions, 
ZPG  and  APG  flowfield  calculations  are  extremely  sensitive  to  such  variations. 

3.  Favorable  pressure  gradients  tend  to  damp  out  turbulent  fluctuations  in  compressible 
boundary  layer  flow.  The  effect  of  the  imposition  of  a  favorable  pressure  gradient  on 
the  flowfields  studied  was  to  reduce  the  turbulent  shear  stress  and  alter  the  distribu¬ 
tion  of  the  turbulent  kinetic  energy  through  the  boundary  layer.  Wall  shear  stress  is 
also  reduced  in  regions  of  favorable  pressure  gradient. 

4.  Adverse  pressure  gradients  tend  to  increase  the  magnitude  of  the  turbulent  fluctu¬ 
ations  in  compressible  boundary  layer  flowfields.  The  effect  of  an  adverse  pressure 
gradient  on  the  flowfields  studied  was  to  increase  the  turbulent  shear  stress  through 
the  boundary  layer.  Additionally,  the  peak  turbulent  kinetic  energy  and  the  wall 
shear  stress  are  increased  under  the  influence  of  the  adverse  pressure  gradient. 
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6.2  Recommendations  for  Further  Investigation 


1.  Decrease  the  solution  time  by  enabling  the  implicit  solution  of  the  k  —  u  equations. 
If  running  in  implicit  mode,  the  flow  solver  in  its  present  form  utilizes  Symmetric 
Gauss-Seidel  Simultaneous  Over-Relaxation  (SOR)  by  Lines  to  rapidly  advance  the 
solution  utilizing  CFL  numbers  as  high  as  1.0x10®.  In  this  mode,  converged  laminar 
solutions  are  possible  in  as  few  as  200  iterations.  Enabling  the  solution  of  the  k 
and  u  equations  to  proceed  at  this  rate  would  significantly  reduce  computation  time, 
dramatically  increasing  this  code’s  utility. 

2.  Move  on  to  more  challenging  flowfields.  Now  that  the  performance  of  this  algo¬ 
rithm/turbulence  model  combination  has  been  documented  for  compressible  turbu¬ 
lent  flows  involving  mild  pressure  gradients,  evaluate  its  applicability  to  separated 
flows  or  flows  over  more  complex  geometries.  Specifically,  the  supersonic  slot  injection 
flowfields  studied  by  Tucker  (23)  stand  out  as  excellent  test  cases  for  an  investigation 
of  this  nature. 

3.  Integrate  additional  sophisticated  turbulence  models.  For  the  mild  pressure  gradi¬ 
ent  cases  presented  in  this  study,  the  Baldwin-Lomax  turbulence  model  is  relatively 
well  behaved  and  helps  to  provide  a  baseline  against  which  the  performance  of  the 
k  —  u  model  may  be  Judged.  For  the  massively  separated  flows  mentioned  above, 
algebraic  models  such  as  Baldwin-Lomax  tend  not  to  fare  very  well.  The  addition 
of  a  Reynolds-Stress  turbulence  model,  directly  solving  the  Reynolds-stress  equa¬ 
tion,  will  provide  additional  capabilities  for  computations  involving  more  challenging 
flowfields.  Failing  that,  the  addition  of  one  or  several  other  two-equation  models  will 
enhance  this  code’s  utility  as  an  analysis  tool  and  serve  as  a  basis  for  compjirison  for 
investigations  involving  these  more  complicated  flowfields. 
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Appendix  A.  Strong  Conservation  Law  form  of  Governing  Differential  Equations 

This  section  presents  the  tramsformation  of  the  governing  equations  from  Cartesian 
form  to  the  Strong  Conservation  Law  form.  The  governing  partial  differential  equations  in 
vector  form  in  Cartesian  coordinates 


dU  dF  dG  „ 

- L  -  -I-  -  =  0 

di  dx  dy 


(A.l) 


may  be  recast  using  the  chain  rule  of  partial  differentiation 


dx  dr, 


(A.2) 


dy  dr, 

into  the  Weak  (or  Chain  Rule)  Conservation  form  of  the  governing  equations 


(A.3) 


du  ^  dF  dF  ^  dG  dG  ^ 
■^  +  =  0 


(A.4) 


The  next  step  in  the  transformation  is  to  divide  the  above  equation  by  the  Jacobian  of  the 
forward  transformation  and  regroup 


J  dt^\J  di^  J  di)^\j  dr,'^  J  dr,) 
The  first  term  in  parentheses  on  the  LHS  may  be  re-expressed  as 


(A.5) 


Similarly,  the  second  term  in  parentheses  may  be  expressed  as 


(A.6) 


7^'  + =  (j^  +  f «)  -  (t).  ^  -  (j). « 
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Inserting  these  two  terms  into  Equation  A.5  and  re-grouping  yields 


(A.8) 

Re-expressing  the  terms  in  square  brackets  in  terms  of  the  metrics  of  the  inverse  transfor¬ 
mation  (see  Appendix  D)  reveals  that 

(A.9) 

(A.IO) 

Finally,  the  strong  conservation  law  form  of  the  governing  equations 

is  obtained 

jat^K  j  +  l  j  J,-" 

(A.ll) 

Now,  defining 

+ 

II 

(A.12) 

and 

II 

(A.13) 

the  strong  conservation  law  form  of  the  governing  equations  may  be  rewritten  as 


^  -h  Ef  -h  G,  =  0  (A.14) 
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Appendix  B.  Finite  Volume  Method  Overview 


This  section  presents  an  overview  of  the  Finite  Volume  method,  based  on  the  presen¬ 
tation  given  in  Buter  (3).  Consider  the  governing  differential  equations  of  fluid  mechanics 
in  strong  conservation  form 


J  dt  ^  ^  dr) 


(B.l) 


If  this  equation  is  integrated  over  a  stationary  control  volume  with  vertices  abed,  as  shown 
in  Figure  B,  this  equation  becomes 


Applying  Green’s  Theorem  (13)  to  the  right  hand  side  of  the  above  equation  yields 

if  (Fdrj-GdO  (B.3) 

"  Jabed  \  Ot  J  Jabed 

where  the  line  integral  is  evaluated  by  integrating  counterclockwise  around  the  cell.  Refer¬ 
ring  to  Figure  B,  Equation  B.3  may  be  represented  by  the  following  finite  approximation 


1  /l/n+l  _ 

(B.4) 

—  [FsApat  -1-  FsAijic  -)-  Ff/Arjed  +  FwApda] 

+  +  GsA^ic  A  Gf/A^cd  + 

where  Fe  refers  to  the  evaluation  of  F  at  point  E  and  so  on.  Referring  once  again  the 
uniform,  unitary  computational  mesh  in  Figure  B,  it  is  clear  that 


Afjab  —  AfJ^d  —  A^ic  —  A^jg  —  0 


Arjic  =  A^gi  =  1 
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Figure  B.l  Finite  Volume  Cell  in  Computational  Domain 


^Vda  —  —  1 


which,  when  inserted  into  Equation  B.4  yields 


1  /[/"+>  -U”\  ^ 

J  ^  J  —  +  Gn  —  Fe  —  Fw  (B-5) 


Finally,  defining  6U  =  !/"+'  -  U”  and  multiplying  through  by  J6t  yields 


6U  —  JSt  [(j5  —  Gff  +  Fw  —  Fe]  (B-fi) 
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Appendix  C.  Extrapolation  and  Limiting 


In  cell  centered,  finite  volume  codes,  no  flowfield  information  is  stored  at  the  cell 
interfaces.  Finite  differences  involving  the  fluxes  F  and  G,  however,  require  the  evaluation 
of  these  fluxes  at  those  interfaces.  Therefore,  approximations  of  the  flowfield  properties  at 
the  interfaces  are  required.  A  first  order  approximation  assumes  that  the  flux  at  the  inter¬ 
face  is  the  same  as  the  flux  at  the  cell  center.  A  second  order  approximation  extrapolates 
interface  values  using  the  value  of  the  flux  vector  at  the  cell  center  and  that  of  the  cell 
immediately  upstream.  These  two  methods  of  extrapolation  are  shown  for  the  left  state 
in  Figure  C.l.  From  this  figure  it  is  obvious  that  second  order  extrapolation  utilizes  more 
flowfield  information  and  therefore  should  provide  better  results.  Unfortunately,  second 
order  extrapolation  breaks  down  at  local  maxima  and  minima,  where  extrapolation  based 
on  upstream  property  variance  can  induce  overshoots  at  cell  interfaces.  Thus,  any  extrap¬ 
olation  greater  than  first  order  must  be  controlled,  or  limited,  to  prevent  oscillations  in 
the  vicinity  of  local  maxima  or  minima.  Yee  (28)  developed  the  MINMOD  limiting  scheme 
to  address  this  difficulty.  MINMOD  is  neither  first  nor  second  order,  but  something  in 
between,  providing  most  of  the  benefits  of  second  order  extrapolation  everywhere  but  at 
local  minima  and  maxima,  where  it  reverts  to  first  order.  This  section  describes  MINMOD 
limiting  as  implemented  by  Gaitonde  (8)  within  FPFV. 


C.l  Left  State 


Figure  C.l  First  and  Second  Order  Extrapolation,  Left  State 


Figure  C.2  Sign  Convention  for  Extrapolation  and  Limiting,  Left  State 


If  0  is  a  flowfield  variable  to  be  extrapolated,  the  left  state  is  computed  using  using 
The  sign  convention  for  extrapolation  of  the  left  state  is  shown  in 
Figure  C.2,  with  A  and  B  positive  for  the  configuration  shown,  or 


B  —  ~  4>iJ 

A  = 


(C.l) 


In  general,  not  including  the  case  where  A  =  B  =  0,  six  cases  must  be  considered  and 
handled  differently,  depending  on  the  signs  and  magnitudes  of  A  and  B.  All  six  cases  are 
handled  within  FPFV  rather  elegantly  by  the  following  functioned  evaluation 


(C.2) 


where  FL1(A,B)  is  a  function  defined  by 


FL1(A,B)  =  sign(l,A)max(0,min(|A|,Bsign(l,A)))  (C.3) 

where 


|A|:  B>0 
-|A|:  B<0 


(C.4) 


These  functional  evaluations  are  perhaps  better  understood  if  presented  in  a  graphical 
manner,  as  is  done  below  for  the  six  cases  under  consideration 


Case  1  :  A  >  B  >  0  :  See  Figure  C.3,  .  =  4>i,i  +  y 

Case  2  :  B  >  A  >  0  :  See  Figure  C.4,  .  =  <f>ij  +  y 
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Figure  C.3  A  >  B  >  0,  Left  State 


CeB  iMerfaces  - 

CeB  Cdittted  State  # 

Extrapolaied  State  ■ 


i-l  i 

Figure  C.4  B  >  A  >  0,  Left  State 

Case  3  :  A,  B  <  0,  |A|  >  |Bl  ;  See  Figure  C.5,  ^ 

Case  4  :  A,  B  <  0,  |B|  >  |A|  :  See  Figure  C,6,  -  4^ 

Case  5  :  A  >  0,  B  <  0  :  Local  Maximum.  See  Figure 

Case  6  :  A  <  0,  B  >  0  :  Local  Minimum.  See  Figure  C.8,  j  =  <l>ij 

C.2  Right  State 

Still  using  <t>  as  the  fiowfield  variable  to  be  extrapolated,  the  right  state  is  computed 
using  <l>ij,  and  The  sign  convention  for  extrapolation  of  the  right  state  is 

shown  in  Figure  C.9,  with  A  and  B  positive  for  the  configuration  shown,  or 

Call  biofacea  — - 

CeOCcaaefed  State  # 

EanpolaiedSaie  ■ 


i-l  i  M 

Figure  C.5  A,  B  <  0,  |A)  >  jB],  Left  State 
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Figure  C.6  A,  B  <  0,  |B|  >  |A|,  Left  State 
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Figure  C.7  Local  Maximum,  Left  State 
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Figure  C.8  Local  Minimum,  Left  State 
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Figure  C.9  Sign  Convention  for  Extrapolation  and  Limiting,  Right  State 
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Cell  Interfaces 
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Figure  C.IO  A  >  B  >  0,  Right  State 


Cell  Intefftcts  — - 
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Figure  C.ll  B  >  A  >  0,  Right  State 


A  —  <t>i+2,j  — 

B  =  -  <i>ij 


(C.5) 


Once  again,  six  cases  must  be  considered  and  handled  differently,  depending  on  the  signs 
and  magnitudes  of  A  and  B.  The  functional  evaluation  presented  in  the  previous  section 
is  used  with  the  following  modification 


-  *t>i+i,}  - 


FL1(A,B) 

2 


(C.6) 


where  FL1(A,B)  is  as  defined  in  Equation  C.3.  Graphically,  the  six  cases  under  consider¬ 
ation  are  as  shown  below 

Case  1  :  A  >  B  >  0  :  See  Figure  C.IO,  -  y 

Case  2  :  B  >  A  >  0  :  See  Figure  C.ll,  =  ^i+ij  -  y 
Case  3  :  A,  B  <  0,  |A|  >  |B|  :  See  Figure  C.12,  =  <i>ij  +  ^ 

Case  4  :  A,  B  <  0,  |B|  >  |A|  :  See  Figure  C.13,  +  *y^ 

Case  5  ;  A  >  0,  B  <  0  :  Local  Maximum.  See  Figure  C.14,  —  4>i+ij 

Case  6  ;  A  <  0,  B  >  0  :  Local  Minimum.  See  Figure  C.15,  =  ^i+ij 
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Cell  laterfaces  — 

Cell  Ce&ieted  Stile  # 

ExtnpoUted  Stue  ■ 


i  i+1  M 

Figure  C.12  A,  B  <  0,  |A|  >  |B|,  Right  State 


Figure  C.13  A,  B  <  0,  |B|  >  |A|,  Right  State 


CeO  bttifaces  <— 

Cell  Cctteied  Swte  • 

Eitnpolaied  Sttte  ■ 


1  M  M 


Figure  C.14  Local  Maximum,  Right  State 


CeDbaofkes  — ~ 

CeS  Cemered  Sate  • 

Extnpolaied  Sate  ■ 


i  Itl  M 


Figure  C.15  Local  Minimum,  Right  State 
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Appendix  D.  Coordinate  Transformation 

Numerically  solving  the  above  system  of  differential  equations  is  most  often  facilitated 
by  transforming  the  frequently  irregular  physical  domain  into  a  rectilinear  computational 
domain  using  evenly  spaced  grid  points.  This  section  presents  the  details  of  a  general 
coordinate  transformation  as  described  by  Buter  (2). 

Thus,  the  physical  {x,y)  space  is  mapped  into  the  computational  space  using 
the  following  coordinate  transformation 

l  =  (D.l) 


V  =  riix,y) 


(D.2) 


Applying  the  chain  rule  of  partial  differentiation  to  the  above  transformation  yields 


^  -c  ^  ^  ^ 

dx 


(D.3) 


(D.4) 


For  the  forward  transformation,  the  total  derivatives  may  be  expressed  as 


d^  =  ^idx  +  ^fdy 


(D.5) 


or 


dp  =  p^dx  +  Pfdy 


(D.6) 


»  - 

dx 

dp 

Px 

Ps 

dy 

(D.7) 
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where  the  Jacobian  of  the  forward  transformation  may  be  defined  as 


(D.8) 


Similarly,  for  the  inverse  transformation,  the  total  derivatives  may  be  expressed  as 


dx  =  X(d^  +  x,,dfj 


(D.9) 


dy  =  y(d^  +  y„dv 


(D.IO) 


If  I, 

yv 


(D.ll) 


Combining  Equations  D.7  and  D.ll  yields 


x„ 


y.?  Vs  Vy 


(D.12) 


If  the  Jacobian  of  the  forward  transformation  is  defined  as 


J  =  ^sVy-  Vxiy 


(D.13) 


the  relationships  between  the  metrics  of  the  forward  and  inverse  transformations  as 


ix  =  Jy„ 


(D.14) 


^y  —  dx^ 


(D.15) 
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=  -JVi 


(D.16) 


r]y  =  Jx^  (D.17) 

Similar  manipulations  reveal  that  the  Jacobian  of  the  inverse  transformation  may  be  de¬ 
fined  as 


J  -  x^Vr,  -  y^Xy 


(D.18) 


and  that  the  Jacobi ans  of  the  forward  and  inverse  transformations  are  related  by 


(D.19) 
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Appendix  E.  Time-Averaged  Navier-Stokes  Equations 

This  appendix  presents  an  overview  of  Reynolds  and  Favre  averaging  as  presented 
by  Bowersox  (1)  and  the  forms  of  the  Navier-Stokes  equations  obtained  through  the  use  of 
each  technique. 

E.l  Reynolds  Averaging 

In  Reynolds  averaging,  flow  variables  are  replaced  by  the  sum  of  a  mean  value  and  a 
fluctuating  component,  as  shown  below  for  the  “random”  flow  variable  (f>. 


=  +  (E.l) 

After  these  substitutions  have  been  made,  the  governing  equations  are  expanded  and  then 
time  averaged.  The  time  average  of  fluctuating  components  is,  by  definition,  zero,  so  terms 
involving  the  product  of  a  mean  and  a  fluctuating  component  disappear,  products  of  two 
or  more  fluctuation  terms  do  not  necessarily  go  to  zero  physically,  so  they  are  preserved 
through  the  averaging  process.  Consider  the  Continuity  equation  in  Cartesian  coordinates 


dp  dpu  dpv 
dt  ^  dx  ^  dy 


(E.2) 


after  substituting  in  the  mean  and  fluctuating  components  for  all  flow  variables  and  ex- 
pajiding,  the  above  equation  becomes 


^{P  +  P')  9  {pu  -I-  pu'  -I-  p'u  A  p'u')  d{pv^-  pv'  -1-  p'v  -t-  p'v') 

at  dx  dy  ■  ° 

Time  averaging  the  resulting  equation  yields  the  Reynolds-averaged  form 


(E.3) 


dp  d  {pu  -b  p'u*)  d  {pv  +  p'v') 
which  may  be  written  in  indicia!  notation  as 


(E.4) 
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(E.5) 


dt  ^  dxj 


Similar  manipulations,  when  performed  on  the  momentum  equations  (22),  yield 


d  {pUj  +  p'u'^  dp^.^.  _  dTjj  +  tT. 

dt  dxj  dxi  dXj 

where  the  Reynolds  stress  tensor  r?  contains  the  additional  fluctuation  terms  added  to  the 
momentum  equations.  In  indicia!  notation,  the  Reynolds  stress  tensor  takes  the  following 
form: 


^Ij  =  -  (/>«<«)  +  “.P'tij  + 


(E.7) 


The  Reynolds  averaged  form  of  the  Energy  equation  (22)  in  terms  of  the  total  enthalpy  is 


d  (pho  +  p'h'^)  Qph^fj.  _  d  (ujTij  +  u't<j  -  qj  -  qj) 
dt  ^  dXj  dxj 


(E.8) 


where  the  turbulent  heat  flux  is  defined  as 


q[  =  +  hop'u'i  +  Uip'h'^  +  p'h'„u'i 


(E.9) 


E.2  Favre  Averaging 

The  form  of  Reynolds  stress  tensor  generated  by  Reynolds  averaging  contains  four 
terms  that  must  be  modeled  if  the  system  of  equations  is  to  be  closed.  In  an  attempt 
to  reduce  the  number  of  terms  that  must  be  modeled,  a  closer  look  at  the  Navier-Stokes 
Equations  suggests  that  mass  averaging  the  flow  variables  is  another  way  to  account  for 
the  turbulent  fluctuations  within  the  flow.  The  premise  of  replacing  flow  quantities  with  a 
mean  plus  a  fluctuating  component  does  not  change,  but  the  rules  by  which  the  averaging 
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process  eliminates  combinations  of  terms  changes  commensurate  with  the  philosophy  of 
mass,  or  Favre  ,  averaging.  Consider  once  .again  a  random  flow  variable,  (f>‘. 

=  (E.IO) 


where 


(E.11) 


In  Favre  averaging,  the  averaging  process  eliminates  only  those  terms  that  involve  the 
product  of  the  mean  density  and  a  fluctuating  term.  This  feature  arises  from  the  following 
set  of  manipulations 


p(f>  =  p<f/' 


(E.12) 


Averaging  yields 


p4>  =  +  p<(/' 

=  p^  +  p^ 

which,  when  combined  with  Equation  E.ll,  yields 


(E.13) 


p<P  =  p!^  +  p(f/> 

P 


(E.14) 


or 


p4>"  =  0  (E.15) 

Applying  these  rules  to  the  governing  differential  equations  (22)  yields  the  Favre 
averaged  form  of  the  Continuity  Equation 
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dp  dp^  ^ 

at  dxi 


(E.16) 


the  Momentum  Equations, 


djpuj)  dpUjUj  _  dTjj  +  tT. 

at  axj  axi  axj 


(E.17) 


where  the  Reynolds  stress  tensor  tT.  is  defined  in  indicia!  notation  as 


'I'lj  =  -/>« 


(E.18) 


and  the  Energy  equation 


d{ph)  dphuj  ^  ap  ap 
at  axj  at  ^  axj 


=  ^  +  u.2t  +  u'!^  4- 

at  ^  axj  axj  ax,  ax^ 


(E.19) 


where  the  turbulent  heat  flux  is  given  by 


qj  =  ph"uy 


(E.20) 


Notice  that  the  Reynolds  stress  tensor  obtained  through  Favre  Averaging  contains  three 
fewer  terms  than  that  obtained  from  Reynolds  averaging.  In  fact,  the  Favre  averaged  term 
looks  much  like  the  “incompressible”  portion  of  the  Reynolds  averaged  term,  as  would  be 
obtained  by  setting  p'  =  0.  The  following  exercise  more  formally  shows  the  relationship 
between  the  Favre  -averaged  and  Reynolds- averaged  terms.  The  instantaneous  velocity 
may  be  represented  in  terms  of  both  Reynolds-averaged  and  Favre  -averaged  variables 


«  =  «  -f  «'  =  « -f  u" 


(E.21) 
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V  =  u  -|-  v'  =  w  +  v" 

Solving  for  the  '  quantities, 

v!  =  u  —  u-\-  u"  (E.22) 

v'  =  V  —  V  v" 

Next,  take  the  time-average  the  Favre  portions  of  Equation  E.21 

u  +  lF  (E.23) 

V  -h  v" 

and  use  them  to  replace  the  u  and  v  terms  in  Equation  E.22 

u'  —  u  u"  —  u  —  u"  —  u"  —  u"  (E.24) 

v'  =  V  v"  —  V  —  v"  =  v"  —  v" 

Next,  take  the  product  of  these  two  terms  and  multiply  the  result  by  the  mean  density 

pu'v'  =  p  v,"v"  —  p  u"  v"  —  p  v"  u"  -f-  p  u"  v"  ~  p  v,"v"  —  p  u"  v"  (E.25) 

Then,  remembering  the  principles  of  Favre  averaging  in  terms  of  a  general  flow  variable  <f> 


u  =  u-\-u"  — 
V  z=  V  v"  = 
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(E.26) 


*f>  =^  +  (f>" 

0"  =  <f>  —  ^ 


the  Favre  fluctuating  component  may  be  re-written  as 


(E.27) 


This  result  may  be  inserted  into  the  right  hand  side  (RHS)  of  Equation  E.25  to  yield 


_  _  p'u'  o'v' 

pu'v'  =  pu"v"  -j — z — ^ 
P  P 


(E.28) 


If  the  second  term  on  the  RHS  of  this  equation  is  negligible,  as  may  be  expected  of  a 
fourth-order  fluctuation  term,  then  the  above  relation  reduces  to 


pu'v'  «  pu"v"  (E.29) 

which  supports  the  observation  made  previously.  Thus,  solving  the  Favre  averaged  form 
of  the  governing  equations  neglects  the  effects  of  p',  manifested  by  the  elimination  of 
the  fourth-order  fluctuation  term  in  the  above  expression.  Justification  for  doing  so  lies 
in  Morkovin’s  Hypothesis  (17),  which  states  that  the  effect  of  density  fluctuations  on 
turbulent  flow  may  be  neglected  as  long  as  the  magnitude  of  the  density  fluctuations  is 
small  compared  to  the  mean  density.  Miller  (16)  found  that  density  fluctuations  in  the 
mean  flow  were  on  the  order  of  1%,  but  within  the  boundary  layer  the  magnitude  of 
these  fluctuations  increased  to  roughly  10%.  This  observation  indicates  that  the  effects  of 
density  fluctuations  are  not  negligible  and  should  therefore  technically  be  considered  for 
this  flowfleld.  Modeling  constraints,  however,  dictate  that  the  Favre  averaged  form  of  the 
governing  equations,  which  by  definition  neglect  density  fluctuations,  be  used. 


E-6 


Appendix  F.  AFIT  Supersonic  Pressure- Vacuum  Wind  Tunnel 


The  coordinates  for  the  contour  of  the  AFIT  Supersonic  Pressure- Vacuum  Wind 
Tunnel,  beginning  in  the  upstream  settling  chamber,  are: 


Ceiling 


X  (in) 

y  (in) 

X  (in) 

y  (in) 

O.OOOOE-i-00 

3.2910E+00 

2.2240E-t-00 

2.2240E-I-00 

1.5000E-f-00 

3.2910E-I-00 

2.2590E+00 

2.1420E-I-00 

1.5340E-f-00 

3.2870E-1-00 

2.2930E-I-00 

2.0600E-f00 

1.5690E-I-00 

3.2770E-i-00 

2.3280E-I-00 

1.9770E-i-00 

1.6030E-I-00 

3.2610E-1-00 

2.3620E-I-00 

1.8940E-f-00 

1.6380E-f00 

3.2390E-I-00 

2.3970E-f00 

1.8110E-f00 

1.6720E-I-00 

3.2120E+00 

2.4310E-I-00 

1.7290E-f00 

1.7070E-f00 

3.1790E-I-00 

2.4660E-b00 

1.6470E-I-00 

1.7410E-I-00 

3.1400E+00 

2.5000E-f00 

1.5670E-}-00 

1.7760E-1-00 

3.0970E-I-00 

2.5350E-I-00 

1.4870E-1-00 

1.8100E-1-00 

3.0500E-I-00 

2.5690E-f00 

1.4090E-t-00 

1.8450E-t-00 

2.9980E-I-00 

2.6040E+00 

1.3330E-I-00 

1.8790E-f00 

2.9420E-f-00 

2.6380E+00 

1.2590E-1-00 

1.9140E-I-00 

2.8830E-I-00 

2.6730E-t-00 

1.1870E-I-00 

1.9480E-I-00 

2.8200E-1-00 

2.7070E-I-00 

1.1180E-f00 

1.9830E-f00 

2.7530E-I-00 

2.7420E-1-00 

1.0520E+00 

2.0170E-I-00 

2.6840E-f00 

2.7760E-I-00 

9.8890E-01 

2.0520E-f00 

2.6120E+00 

2.8110E-I-00 

9.2930E-01 

2.0860E-I-00 

2.5380E-I-00 

2.8450E-i-00 

8.7350E-01 

2.1210E-I-00 

2.4620E-f00 

2.8800E-f00 

8.2160E-01 

2.1550E-I-00 

2.3840E-f-00 

2.9140E-I-00 

7.7410E-01 

2.1900E-I-00 

2.3050E-i-00 

2.9490E-1-00 

7.3110E-01 

Ceiling  (continued) 


X  (in) 

y  (>n) 

x(in) 

y  (in) 

2.9830E+00 

6.9290E-01 

5.1000E+00 

1.2600E+00 

3.0180E+00 

6.5990E-01 

5.3030E+00 

1.3310E+00 

3.0520E+00 

6.3220E-01 

5.5120E+00 

1.4010E+00 

3.0870E+00 

6.1020E-01 

5.7290E+00 

L4700E+00 

3.1210E+00 

5.9420E-01 

5.9560E+00 

1.5390E+00 

3.1560E+00 

5.8430E-01 

6.1940E+00 

1.6070E+00 

3.1900E+00 

5.8100E-01 

6.4440E+00 

1.6740E+00 

3.2620E+00 

5.8360E-01 

6.7070E+00 

1.7410E+00 

3.3050E+00 

5.8680E-01 

6.9840E+00 

1.8070E+00 

3.3490E+00 

5.9120E-01 

7.2770E+00 

1.8720E+00 

3.3920E+00 

5.9700E-01 

7.5880E+00 

1.9360E+00 

3.4360E+00 

6.0410E-01 

7.9170E+00 

2.0000E+00 

3.4810E+00 

6.1260E-01 

8.2660E+00 

2.0620E+00 

3.5260E+00 

6.2260E-01 

8.6370E+00 

2.1210E+00 

3.5710E+00 

6.3420E-01 

9.0310E+00 

2.1780E+00 

3.6180E+00 

6.4740E-01 

9.4510E+00 

2.2340E+00 

3.6650E+00 

6.6240E-01 

9.8980E+00 

2.2850E+00 

3.7130E+00 

6.7920E-01 

1.0370E+01 

2.3330E+00 

3.7620E+00 

6.9800E-01 

1.0880E+01 

2.3770E+00 

3.8120E+00 

7.1900E-01 

1.1430E+01 

2.4160E+00 

3.8380E+00 

7.3040E-01 

1.2010E+01 

2.4480E+00 

4.0860E+00 

8.4690E-01 

1.2630E+01 

2.4740E+00 

4.5070E+00 

1.0320E+00 

1.3290E+01 

2.4920E+00 

4.7050E+00 

1.1120E+00 

1.4000E+01 

2.5000E+00 

4.9020E+00 

1.1870E+00 

Floor 


X  (in) 

y  (in) 

x(in) 

y  (in) 

.OOOOE+00 

-.2710E+01 

.2397E+01 

-.1230E+01 

.1500E+01 

-.2710E+01 

.2431E+01 

-.1148E+01 

.1534E+01 

-.2706E+01 

.2466E+01 

-.1066E+01 

.1569E+01 

-.2696E+01 

.2500E+01 

-.9857E+00 

.1603E+01 

-.2680E+01 

.2535E+01 

-.9062E+00 

.1638E+01 

-.2658E+01 

.2569E+01 

-.8283E+00 

.1672E+01 

-.2631E+01 

.2604E+01 

-.7521E+00 

.1707E+01 

-.2598E+01 

.2638E+01 

-.6780E+00 

.1741E+01 

-.2559E+01 

.2673E+01 

-.6063E+00 

.1776E+01 

-.2516E+01 

.2707E+01 

-.5372E+00 

.1810E+01 

-.2469E+01 

.2742E+01 

-.4709E+00 

.1845E+01 

-.2417E+01 

.2776E+01 

-.4079E+00 

.1879E+01 

-.2361E+01 

.2811E+01 

-.3483E+00 

.1914E+01 

-.2302E+01 

.2845E+01 

-.2925E+00 

.1948E+01 

-.2239E+01 

.2880E+01 

-.2406E+00 

.1983E+01 

-.2172E+01 

.2914E+01 

-.1931E+00 

.2017E+01 

-.2103E+01 

.2949E+01 

-.1501E+00 

.2052E+01 

-.2031E+01 

.2983E+01 

-.1119E+00 

.2086E+01 

-.1957E+01 

.3018E+01 

-.7888E-01 

.2121E+01 

-.1881E+01 

.3052E+01 

-.5122E-01 

.2155E+01 

-.1803E+01 

.3087E+01 

-.2923E-01 

.2190E+01 

-.1724E+01 

.3121E+01 

-.1317E-01 

.2224E+01 

-.1643E+01 

.3156E+01 

-.3339E-02 

.2259E+01 

-.1561E+01 

.3190E+01 

.OOOOE+00 

.2293E+01 

-.1479E+01 

.3290E+01 

•OOOOE+OO 

.2328E+01 

,2362E+01 

-.1396E+01 

-.1313E+01 

.1400E+02 

.OOOOE+00 

s 


0  , 


10* 


10*  10* 


Figure  G.8 
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u+  (van  Driest)  vs  y+,  Favorable  Pressure  Gradient  Section,  FS  71.5,  Stream- 
wise  Node  Count  Variation 


71.5,  Streamwise  Node  Count  Variation 
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Figure  G.ll  NondimensionaJ  Turbulent  Shear  Stress  Profile,  Upstream  Flat  Plate  Test 
Section,  FS  44.0,  Streamwise  Node  Count  Variation 


Figure  G.12  NondimensionaJ  Turbulent  Shear  Stress  Profile,  Favorable  Pressure  Gradient 
Section,  FS  71.5,  Streamwise  Node  Count  Variation 
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71.5,  Wall  Spacing  Variation 
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